Wat is het verband tussen vertragingen en het weer?

Dit notebook is gemaakt voor het AI4 semester van Fontys ICT, als onderdeel van een eigen uitdaging.

Het is geschreven door mij, Anne Koppers.

In dit notebook ga ik proberen te onderzoeken of het weer invloed heeft op vertragingen in het openbaar treinverkeer.

Met behulp van deze datasets hoop ik te kunnen ontdekken of en waar er een sterke correlatie is.

In mijn uitdaging ga ik gebruik maken van Python Pandas, een bibliotheek om datasets in te laden. Daarnaast ga ik ook gebruik maken van Plotly om de data op een visueel aantrekkelijke manier weer te geven, en Scikit-learn (sklearn) voor de machine learning algoritmes.

Afhankelijkheden voor de uitdaging:

# dependances
import os
import numpy as np
# import matplotlib as mpl
from matplotlib import pyplot as plt
import seaborn as sns
import pandas as pd
from datetime import datetime
# from geopy.exc import GeocoderTimedOut
# from geopy.geocoders import Nominatim
# from geopy.distance import geodesic
import matplotlib.image as mpimg

import warnings
# warnings.filterwarnings('ignore')

from sklearn.impute import KNNImputer

datafolder = 'data/'

max_data = datetime(2023, 4, 1, 0, 0, 0)
calculate_large_data = False

📦 Data opzet

Ik ben op zoek gegaan naar datasets die relevant kunnen zijn voor mijn opdracht.

Dit zijn de datasets die ik heb gekozen en waarom ik deze heb gekozen:

• rijden de treinen. Treinarchief
  deze dataset heeft alle data van het trein verkeer. waaronder de verwachte aankomst tijd en hoeveel minuten hij te laat was.

• rijden de treinen. treinstoringen
  deze dataset heeft informatie over wat voor storingen er vermeld werden. deze will ik gebruiken om te controleren wat mijn vindingen zijn aan het einde.

• rijden de treinen. treinstations
  Deze dataset bevat de coördinaten van de trein stations van de dataset.
  de datasets van Rijden de treinen woord gebruikt door nieuws programmas and kranten zoals AD en Volkskrant. En heeft data die ik nodig heb om mijn berekeningen te doen Dit maakt ze betrouwbaar en geschrikt voor wat ik will berijken.

• KNMI weerstations Uurwaarnemingen
  Dit is een dataset van de overheid. dit zoorgt ervoor dat het betrouwbaar genoeg is om er van overtuigt the zijn dat het betrouwbaar is. dus ik kan deze dataset gebruiken om te bepalen wat het weer was dichtbij een trein station.

In deze datasets is informatie te vinden over al het treinverkeer, het treinverkeer dat als storing is aangemerkt, het weer per uur en de locatie van het treinstation.

🚄 Laden van het Treinarchief

De onderstaande code laadt het treinarchief in. Volgens Rijden de Treinen heeft dit treinarchief de volgende kolommen:

• Service:RDT-ID Dit ID wordt gebruikt voor deze rit. Het heeft geen nuttige betekenis buiten het uniek identificeren van een enkele rit op een enkele datum.

• Service:Date De geplande datum voor deze rit.

• Service:Type Het soort trein, zoals Intercity, Sprinter of ICE International.
• Service:Company Vervoerder die deze rit uitvoert, zoals NS of Arriva.
• Service:Train number Het treinnummer (rit-nummer) voor deze rit identificeert deze rit uniek op deze datum. Dit nummer wordt soms ook gecommuniceerd naar reizigers (vooral voor internationale treinen).

Een enkele rit kan soms meerdere treinnummers hebben. Bijvoorbeeld wanneer een trein gesplitst wordt in twee delen, of wanneer een trein op een groot station halverwege van treinnummer verandert.

De rit-datum is niet altijd de logische datum. Bijvoorbeeld: een rit die vertrekt om 23:59 op 31 juli en aankomt om 02:00 op 1 augustus heeft een rit-datum van 31 juli. Vertragingen hebben geen invloed op de rit-datum.

Dit ID zal meerdere keren voorkomen in de CSV-bestanden omdat dit ID uniek is voor iedere rit. Er is ook een kolom Stop:RDT-ID met een uniek ID voor iedere stop.

• Service:Completely cancelled Deze kolom is true wanneer alle stops van deze rit opgeheven zijn. Of met andere woorden: wanneer de trein helemaal niet rijdt.
• Service:Partly cancelled Deze kolom is true wanneer een of meer stops van deze rit opgeheven zijn. Of met andere woorden: wanneer de trein niet rijdt op een deel van het traject.
• Service:Maximum delay De hoogste vertraging (in minuten) van alle stops van deze rit.
• Stop:RDT-ID Unieke ID voor deze stop. Dit ID is uniek voor iedere stop in de dataset. Het heeft verder geen nuttige betekenis.
• Stop:Station code Code (verkorting) van de stationsnaam.
• Stop:Station name De naam van het station.
• Stop:Arrival time Geplande aankomsttijd in RFC 3339-formaat. Deze kolom is leeg wanneer er geen aankomst gepland was.
• Stop:Arrival delay Vertraging bij aankomst in minuten. Deze kolom is leeg wanneer er geen aankomst gepland was.
• Stop:Arrival cancelled Deze kolom is true wanneer de aankomst op deze stop opgeheven is. Deze kolom is leeg wanneer er geen aankomst gepland was.
• Stop:Departure time Geplande vertrektijd in RFC 3339-formaat. Deze kolom is leeg wanneer er geen vertrek gepland was.
• Stop:Departure delay Vertraging bij vertrek in minuten. Deze kolom is leeg wanneer er geen vertrek gepland was.
• Stop:Departure cancelled Deze kolom is true wanneer het vertrek op deze stop opgeheven is. Deze kolom is leeg wanneer er geen vertrek gepland was.

In de code laad ik de dataset in en verander ik de datum naar een datetime-formaat, zodat ik deze kan beperken tot het zoekgebied dat ik momenteel heb in de weerdataset. Vervolgens zet ik 'type' en 'bedrijf' om naar kleine letters, omdat de dataset soms handmatig is aangepast en dan worden hoofdletters vergeten.

Ik gebruik deze dataset om te kijken welke treinen hebben gereden en welke treinen vertraging hadden.

Train_trafic_2023 = pd.read_csv(datafolder + 'services-2023/services-2023.csv', sep=',', low_memory=False)

Train_trafic_2023['Service:Date'] = pd.to_datetime(Train_trafic_2023['Service:Date'], format="%Y-%m-%d")
Train_trafic_2023 = Train_trafic_2023[Train_trafic_2023['Service:Date'] < max_data]

Train_trafic_2023['Service:Type'] = Train_trafic_2023['Service:Type'].str.lower()
Train_trafic_2023['Service:Company'] = Train_trafic_2023['Service:Company'].str.lower()
Train_trafic_2023['Stop:Station name'] = Train_trafic_2023['Stop:Station name'].str.lower()

Train_trafic_2023['Stop:Arrival time'] = pd.to_datetime(Train_trafic_2023['Stop:Arrival time'])
Train_trafic_2023['Stop:Departure time'] = pd.to_datetime(Train_trafic_2023['Stop:Departure time'])

# Train_trafic_2023 = Train_trafic_2023[Train_trafic_2023['Service:Company'] == 'ns']

display(Train_trafic_2023.head())
display(Train_trafic_2023.info())

	Service:RDT-ID	Service:Date	Service:Type	Service:Company	Service:Train number	Service:Completely cancelled	Service:Partly cancelled	Service:Maximum delay	Stop:RDT-ID	Stop:Station code	Stop:Station name	Stop:Arrival time	Stop:Arrival delay	Stop:Arrival cancelled	Stop:Departure time	Stop:Departure delay	Stop:Departure cancelled
0	10367702	2023-01-01	intercity	ns	1410	False	False	3	93067918	RTD	rotterdam centraal	NaT	NaN	NaN	2023-01-01 02:00:00+01:00	1.0	False
1	10367702	2023-01-01	intercity	ns	1410	False	False	0	93067919	DT	delft	2023-01-01 02:13:00+01:00	0.0	False	2023-01-01 02:13:00+01:00	0.0	False
2	10367702	2023-01-01	intercity	ns	1410	False	False	0	93067920	GV	den haag hs	2023-01-01 02:21:00+01:00	0.0	False	2023-01-01 02:22:00+01:00	0.0	False
3	10367702	2023-01-01	intercity	ns	1410	False	False	0	93067921	LEDN	leiden centraal	2023-01-01 02:35:00+01:00	0.0	False	2023-01-01 02:40:00+01:00	0.0	False
4	10367702	2023-01-01	intercity	ns	1410	False	False	0	93067922	SHL	schiphol airport	2023-01-01 02:56:00+01:00	0.0	False	2023-01-01 03:03:00+01:00	0.0	False

<class 'pandas.core.frame.DataFrame'>
Index: 5142250 entries, 0 to 11473990
Data columns (total 17 columns):
 #   Column                        Dtype         
---  ------                        -----         
 0   Service:RDT-ID                int64         
 1   Service:Date                  datetime64[ns]
 2   Service:Type                  object        
 3   Service:Company               object        
 4   Service:Train number          int64         
 5   Service:Completely cancelled  bool          
 6   Service:Partly cancelled      bool          
 7   Service:Maximum delay         int64         
 8   Stop:RDT-ID                   int64         
 9   Stop:Station code             object        
 10  Stop:Station name             object        
 11  Stop:Arrival time             object        
 12  Stop:Arrival delay            float64       
 13  Stop:Arrival cancelled        object        
 14  Stop:Departure time           object        
 15  Stop:Departure delay          float64       
 16  Stop:Departure cancelled      object        
dtypes: bool(2), datetime64[ns](1), float64(2), int64(4), object(8)
memory usage: 637.5+ MB

None

🌥️ Laden van KNMI Uurwaarnemingen

In de onderstaande code laad ik de Uurwaarnemingen van het KNMI in. Dit is de reden waarom het merendeel van het project alleen over de eerste twee maanden van het jaar gaat.

In de weergegevens dataset hebben we:

• STN dit is het weer station nummer
• YYYYMMDD datum (YYYY=jaar;MM=maand;DD=dag)
• start_time tijd (HH=uur; UT.12 UT=13 MET; 14 MEZT. Uurvak 05 loopt van 04.00 UT tot 5.00 UT)
• DD Windrichting (in graden) gemiddeld over de laatste 10 minuten van het afgelopen uur (360=noord; 90=oost; 180=zuid; 270=west; 0=windstil 990=veranderlijk.
• FH Uurgemiddelde windsnelheid (in 0.1 m/s).
• FF Windsnelheid (in 0.1 m/s) gemiddeld over de laatste 10 minuten van het afgelopen uur.
• FX Hoogste windstoot (in 0.1 m/s) over het afgelopen uurvak.
• temp Temperatuur (in 0.1 graden Celsius) op 1.50 m hoogte tijdens de waarneming.
• T10N Minimumtemperatuur (in 0.1 graden Celsius) op 10 cm hoogte in de afgelopen 6 uur.
• TD Dauwpunts-temperatuur (in 0.1 graden Celsius) op 1.50 m hoogte tijdens de waarneming.
• SQ Duur van de zonneschijn (in 0.1 uren) per uurvak; berekend uit globale straling (-1 for <0.05 uur).
• Q Globale straling (in J/cm2) per uurvak.
• DR Duur van de neerslag (in 0.1 uur) per uurvak.
• RH Uursom van de neerslag (in 0.1 mm) (-1 voor <0.05 mm).
• P Luchtdruk (in 0.1 hPa) herleid naar zeeniveau; tijdens de waarneming.
• VV Horizontaal zicht tijdens de waarneming (0=minder dan 100m; 1=100-200m; 2=200-300m;...; 49=4900-5000m; 50=5-6km; 56=6-7km; 57=7-8km; ...; 79=29-30km; 80=30-35km; 81=35-40km;...; 89=meer dan 70km).
• N Bewolking (bedekkingsgraad van de bovenlucht in achtsten); tijdens de waarneming (9=bovenlucht onzichtbaar).
• U Relatieve vochtigheid (in procenten) op 1.50 m hoogte tijdens de waarneming.
• WW Weercode (00-99); visueel(WW) of automatisch(WaWa) waargenomen; voor het actuele weer of het weer in het afgelopen uur.
• IX Weercode indicator voor de wijze van waarnemen op een bemand of automatisch station (1=bemand gebruikmakend van code uit visuele waarnemingen; 2;3=bemand en weggelaten (geen belangrijk weersverschijnsel; geen gegevens); 4=automatisch en opgenomen (gebruikmakend van code uit visuele waarnemingen); 5;6=automatisch en weggelaten (geen belangrijk weersverschijnsel; geen gegevens); 7=automatisch gebruikmakend van code uit automatische waarnemingen).
• M Mist 0=niet voorgekomen; 1=wel voorgekomen in het voorgaande uur en/of tijdens de waarneming.
• R Regen 0=niet voorgekomen; 1=wel voorgekomen in het voorgaande uur en/of tijdens de waarneming.
• S Sneeuw 0=niet voorgekomen; 1=wel voorgekomen in het voorgaande uur en/of tijdens de waarneming.
• O Onweer 0=niet voorgekomen; 1=wel voorgekomen in het voorgaande uur en/of tijdens de waarneming.
• Y IJsvorming 0=niet voorgekomen; 1=wel voorgekomen in het voorgaande uur en/of tijdens de waarneming.

Hier laad ik de weergegevens in en zet ik start_time om naar een datetime. Daarbij combineer ik twee kolommen om een datum en tijd in één kolom te krijgen.

Ik gebruik deze dataset om het weer op een bepaalde datum en bij een specifiek stationsnummer te krijgen.

weather_2023 = pd.read_csv(datafolder + '2023-jan-mar.csv', sep=',', low_memory=False, encoding='utf-8', skiprows=80)

weather_2023['start_time'] = weather_2023['start_time'].astype(str) + ' ' + (weather_2023['HH'] - 1).astype(str) + ':00:00'

weather_2023['start_time'] = pd.to_datetime(weather_2023['start_time'], format="%Y%m%d %H:%M:%S").dt.tz_localize('Europe/Amsterdam', nonexistent='shift_forward')

#  drop columns that are not needed
# weather_2023 = weather_2023.drop(columns=['HH', 'WW', 'IX', 'FH', 'FX', 'TD', 'SQ', 'RH'])
weather_2023 = weather_2023.drop(columns=['HH'])

display(weather_2023.head())
display(weather_2023.info())    

	STN	start_time	DD	FH	FF	FX	temp	T10N	TD	SQ	...	VV	N	U	WW	IX	M	R	S	O	Y
0	209	2023-01-01 00:00:00+01:00	220.0	140.0	150.0	170.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
1	209	2023-01-01 01:00:00+01:00	220.0	130.0	130.0	170.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
2	209	2023-01-01 02:00:00+01:00	240.0	140.0	130.0	200.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
3	209	2023-01-01 03:00:00+01:00	220.0	120.0	130.0	170.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
4	209	2023-01-01 04:00:00+01:00	230.0	130.0	140.0	170.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN

5 rows × 24 columns

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 99360 entries, 0 to 99359
Data columns (total 24 columns):
 #   Column      Non-Null Count  Dtype                           
---  ------      --------------  -----                           
 0   STN         99360 non-null  int64                           
 1   start_time  99360 non-null  datetime64[ns, Europe/Amsterdam]
 2   DD          96861 non-null  float64                         
 3   FH          96860 non-null  float64                         
 4   FF          96861 non-null  float64                         
 5   FX          96860 non-null  float64                         
 6   temp        73440 non-null  float64                         
 7   T10N        12240 non-null  float64                         
 8   TD          73440 non-null  float64                         
 9   SQ          69120 non-null  float64                         
 10  Q           69120 non-null  float64                         
 11  DR          71279 non-null  float64                         
 12  RH          71280 non-null  float64                         
 13  P           51840 non-null  float64                         
 14  VV          53999 non-null  float64                         
 15  N           43200 non-null  float64                         
 16  U           73440 non-null  float64                         
 17  WW          26196 non-null  float64                         
 18  IX          99360 non-null  int64                           
 19  M           53945 non-null  float64                         
 20  R           53945 non-null  float64                         
 21  S           53945 non-null  float64                         
 22  O           53945 non-null  float64                         
 23  Y           53945 non-null  float64                         
dtypes: datetime64[ns, Europe/Amsterdam](1), float64(21), int64(2)
memory usage: 18.2 MB

None

De weer dataset bestaat uit twee datasets. De eerste bevat de gegevens van de weerstations en de tweede gaat over de locaties van de weerstations.

De dataset van de weerstationlocaties heeft de volgende kolommen:

• STN: Dit is het nummer van het weerstation.
• LON(east): Dit is de lengtegraad van het weerstation.
• LAT(north): Dit is de breedtegraad van het weerstation.
• ALT(m): Dit is de hoogte waarop het station zich bevindt.
• NAME: Dit is de naam van de locatie waar het weerstation zich bevindt.

Met deze gegevens hoop ik te kunnen achterhalen wat het weer was op het aankomststation.

De weergegevens zijn gekoppeld aan de locatie op basis van het stationnummer.

weather_station_location = pd.read_csv(datafolder + 'Weather-stations-loc.csv', sep=',', low_memory=False, encoding='utf-8')
display(weather_station_location.head())
display(weather_station_location.info())

	STN	LON(east)	LAT(north)	ALT(m)	NAME
0	209	4.518	52.465	0.0	IJmond
1	210	4.430	52.171	-0.2	ValkenburgZh
2	215	4.437	52.141	-1.1	Voorschoten
3	225	4.555	52.463	4.4	IJmuiden
4	235	4.781	52.928	1.2	DeKooy

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 50 entries, 0 to 49
Data columns (total 5 columns):
 #   Column      Non-Null Count  Dtype  
---  ------      --------------  -----  
 0   STN         50 non-null     int64  
 1   LON(east)   50 non-null     float64
 2   LAT(north)  50 non-null     float64
 3   ALT(m)      50 non-null     float64
 4   NAME        50 non-null     object 
dtypes: float64(3), int64(1), object(1)
memory usage: 2.1+ KB

None

🏟️ Treinstations

De onderstaande dataset bevat alle treinstation Coördinaten.

Ik zal deze dataset gebruiken om de locaties van de treinstations te vinden. Aangezien ik alleen de locaties nodig heb, zal ik alleen die onthouden, naast de belangrijke gegevens zoals mogelijke namen.

Ik hoop deze dataset te kunnen gebruiken om te vergelijken waar een station zich bevindt.

In deze dataset staan de volgende gegevens:

• rdt_id: Intern ID dat door Rijden de Treinen wordt gebruikt.
• code: De unieke afkorting van het station. Meestal is dit een afkorting van de stationsnaam, maar de naam is niet altijd herkenbaar (bijvoorbeeld de meeste stations in Den Haag). Een stationcode is maximaal 8 tekens lang.
• uic: De UIC-code van een station. Dit is een internationaal uniek identificatienummer voor stations. Nederlandse stations hebben een UIC-code die begint met 84 (bijvoorbeeld 8400058 voor Amsterdam Centraal).
• name_short: Korte versie van de stationsnaam. Vaak wordt een groot deel van de stationsnaam weggelaten (bijvoorbeeld Sciencepark in plaats van Amsterdam Science Park).
• name_medium: De middellange versie van de stationsnaam. Deze versie van de stationsnaam wordt vaak gebruikt wanneer er geen ruimte is voor de volledige naam, maar kan ook niet alle naamelementen bevatten (bijvoorbeeld Science Park in plaats van Amsterdam Science Park).
• name_long: De lange (of volledige) stationsnaam.
• slug: Unieke naam voor het station, in kleine letters en zonder spaties. Dit is hetzelfde als de URL van de stationspagina op de Rijden de Treinen website.
• country: Het land waarin het station zich bevindt. NL is Nederland, D is Duitsland, B is België, F is Frankrijk, A is Oostenrijk, CH is Zwitserland, GB is het Verenigd Koninkrijk.
• type: Het soort station. Mogelijke waarden zijn (van groot naar klein station): megastation, knooppuntIntercitystation, intercitystation, knooppuntSneltreinstation, sneltreinstation, knooppuntStoptreinstation, stoptreinstation, facultatiefStation. De categorie knooppunt impliceert dat deze stations bedoeld zijn voor het overstappen op een andere trein. Een facultatiefStation is een station waar normaal gesproken geen treinen stoppen, behalve bij speciale gelegenheden (zoals voetbalwedstrijden).
• geo_lat: Breedtegraad van het station, in decimale graden.
• geo_lng: Lengtegraad van het station, in decimale graden.

train_station_locations = pd.read_csv(datafolder + 'stations-2023-09.csv', low_memory=False, encoding='utf-8')

display(train_station_locations.head())
display(train_station_locations.info())

	id	code	uic	name_short	name_medium	name_long	slug	country	type	geo_lat	geo_lng
0	266	HT	8400319	Den Bosch	's-Hertogenbosch	's-Hertogenbosch	s-hertogenbosch	NL	knooppuntIntercitystation	51.690480	5.293620
1	269	HTO	8400320	Dn Bosch O	's-Hertogenb. O.	's-Hertogenbosch Oost	s-hertogenbosch-oost	NL	stoptreinstation	51.700554	5.318333
2	227	HDE	8400388	't Harde	't Harde	't Harde	t-harde	NL	stoptreinstation	52.409168	5.893611
3	8	AHBF	8015345	Aachen	Aachen Hbf	Aachen Hbf	aachen-hbf	D	knooppuntIntercitystation	50.767800	6.091499
4	818	AW	8015199	Aachen W	Aachen West	Aachen West	aachen-west	D	stoptreinstation	50.780360	6.070715

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 591 entries, 0 to 590
Data columns (total 11 columns):
 #   Column       Non-Null Count  Dtype  
---  ------       --------------  -----  
 0   id           591 non-null    int64  
 1   code         590 non-null    object 
 2   uic          591 non-null    int64  
 3   name_short   591 non-null    object 
 4   name_medium  591 non-null    object 
 5   name_long    591 non-null    object 
 6   slug         591 non-null    object 
 7   country      591 non-null    object 
 8   type         591 non-null    object 
 9   geo_lat      591 non-null    float64
 10  geo_lng      591 non-null    float64
dtypes: float64(2), int64(2), object(7)
memory usage: 50.9+ KB

None

📃 Data voorbeelden

Om een indruk te krijgen van de data, ga ik hieronder 20 items uit elke dataset laden. Dit is om te controleren of de data correct is ingeladen en of de data er bruikbaar uitziet.

🚄 Treinarchief laden

Hier laad ik voorbeelden van het treinverkeer. Daarnaast laat ik ook zien of er NULL-waarden in de dataset zitten. Bij dit treinarchief zijn NULL-waarden niet problematisch, omdat ze een belangrijke betekenis hebben. Een NULL-waarde betekent namelijk dat het betreffende item niet voorkwam. Dus als er bij een vertraging een NULL-waarde staat, betekent dit dat er geen vertraging is gemeld.

sample_size = 20

display(Train_trafic_2023.sample(sample_size))
display(Train_trafic_2023.isna().sum())

	Service:RDT-ID	Service:Date	Service:Type	Service:Company	Service:Train number	Service:Completely cancelled	Service:Partly cancelled	Service:Maximum delay	Stop:RDT-ID	Stop:Station code	Stop:Station name	Stop:Arrival time	Stop:Arrival delay	Stop:Arrival cancelled	Stop:Departure time	Stop:Departure delay	Stop:Departure cancelled
1335763	10506748	2023-01-25	sprinter	ns	5620	False	False	1	94408544	ZL	zwolle	NaT	NaN	NaN	2023-01-25 06:53:00+01:00	0.0	False
4110892	10808118	2023-03-14	sprinter	ns	6651	False	False	0	97183644	O	oss	2023-03-14 14:44:00+01:00	0.0	False	2023-03-14 14:44:00+01:00	0.0	False
1365242	10509707	2023-01-25	sprinter	ns	6654	False	False	0	94438037	HTO	's-hertogenbosch oost	2023-01-25 16:36:00+01:00	7.0	False	2023-01-25 16:36:00+01:00	7.0	False
2014437	10579810	2023-02-04	stopbus i.p.v. trein	ns	981180	False	False	0	95086757	RLB	rotterdam lombardijen	2023-02-04 10:45:00+01:00	0.0	False	NaT	NaN	NaN
2093257	10588656	2023-02-07	stoptrein	r-net	7226	False	False	0	95166072	LDM	leerdam	2023-02-07 09:25:00+01:00	0.0	False	2023-02-07 09:26:00+01:00	3.0	False
4220042	10819178	2023-03-16	sprinter	ns	4841	False	True	0	97292819	ASS	amsterdam sloterdijk	2023-03-16 11:01:00+01:00	0.0	True	2023-03-16 11:01:00+01:00	0.0	True
3082746	10694166	2023-02-24	intercity	ns	13929	False	False	0	96155555	RM	roermond	2023-02-24 10:39:00+01:00	0.0	False	2023-02-24 10:40:00+01:00	0.0	False
283397	10396990	2023-01-06	intercity	ns	3038	False	False	0	93356201	HDRZ	den helder zuid	2023-01-06 14:22:00+01:00	0.0	False	2023-01-06 14:22:00+01:00	0.0	False
325687	10401343	2023-01-07	stoptrein	nmbs	5378	False	False	0	93398491	LUIK	liège-guillemins	2023-01-07 08:52:00+01:00	0.0	False	NaT	NaN	NaN
4116115	10808647	2023-03-14	intercity	ns	8859	False	False	0	97188868	BDG	bodegraven	2023-03-14 16:15:00+01:00	0.0	False	2023-03-14 16:16:00+01:00	0.0	False
822027	10453043	2023-01-16	sprinter	ns	5828	False	False	0	93894828	ASDM	amsterdam muiderpoort	2023-01-16 10:13:00+01:00	1.0	False	2023-01-16 10:13:00+01:00	1.0	False
4056524	10802621	2023-03-13	sprinter	ns	4850	False	True	0	97129289	ASS	amsterdam sloterdijk	2023-03-13 17:30:00+01:00	0.0	True	2023-03-13 17:30:00+01:00	0.0	True
936053	10464615	2023-01-18	sprinter	ns	6329	False	False	0	94008829	GVC	den haag centraal	2023-01-18 08:07:00+01:00	0.0	False	NaT	NaN	NaN
1757546	10552242	2023-02-01	sprinter	ns	5130	False	False	0	94830339	RTD	rotterdam centraal	2023-02-01 10:34:00+01:00	0.0	False	2023-02-01 10:35:00+01:00	0.0	False
654573	10434627	2023-01-13	stoptrein	arriva	37209	False	False	0	93727372	DEI	deinum	2023-01-13 06:19:00+01:00	0.0	False	2023-01-13 06:19:00+01:00	0.0	False
4829088	10883708	2023-03-27	stoptrein	r-net	7214	False	False	3	97901800	GDM	geldermalsen	NaT	NaN	NaN	2023-03-27 06:14:00+02:00	0.0	False
756988	10445510	2023-01-14	snelbus i.p.v. trein	ns	980516	False	False	0	93829547	RSN	rijssen	2023-01-14 23:37:00+01:00	0.0	False	2023-01-14 18:10:00+01:00	0.0	False
4583630	10856164	2023-03-22	sprinter	ns	3364	False	False	0	97656435	HFD	hoofddorp	2023-03-22 19:16:00+01:00	0.0	False	2023-03-22 19:17:00+01:00	0.0	False
2259259	10605578	2023-02-09	sprinter	ns	6986	False	False	0	95332056	CL	culemborg	2023-02-10 00:17:00+01:00	1.0	False	2023-02-10 00:17:00+01:00	1.0	False
4155590	10812674	2023-03-15	sprinter	r-net	8724	False	False	0	97228364	GD	gouda	NaT	NaN	NaN	2023-03-15 09:43:00+01:00	0.0	False

Service:RDT-ID                       0
Service:Date                         0
Service:Type                         0
Service:Company                      0
Service:Train number                 0
Service:Completely cancelled         0
Service:Partly cancelled             0
Service:Maximum delay                0
Stop:RDT-ID                          0
Stop:Station code                 6632
Stop:Station name                    0
Stop:Arrival time               551049
Stop:Arrival delay              551049
Stop:Arrival cancelled          551049
Stop:Departure time             550943
Stop:Departure delay            550943
Stop:Departure cancelled        550943
dtype: int64

🌥️ Laden van KNMI Uurwaarnemingen

Hieronder laad ik voorbeelden van de weerstations en de locaties van de weerstations. Dit heb ik later nodig om te bepalen welk weerstation het dichtst bij het treinstation ligt.

Daarnaast controleer ik hoeveel NULL-waarden er zijn. Dit betekent bij bijna alle kolommen dat een weerstation geen sensor heeft voor die specifieke meting.

T10N is een uitzondering, omdat deze elke 6 uur wordt opgenomen. Daardoor zijn alle andere rijen die ertussen zitten NULL. Dit geldt natuurlijk alleen als het weerstation een temperatuursensor heeft.

display(weather_2023.sample(sample_size))
display(weather_station_location.sample(sample_size))


display(weather_2023.isna().sum())
display(weather_station_location.isna().sum())

	STN	start_time	DD	FH	FF	FX	temp	T10N	TD	SQ	...	VV	N	U	WW	IX	M	R	S	O	Y
33323	273	2023-02-08 11:00:00+01:00	160.0	40.0	50.0	70.0	41.0	-43.0	-33.0	10.0	...	81.0	NaN	58.0	NaN	5	0.0	0.0	0.0	0.0	0.0
84064	348	2023-03-24 16:00:00+01:00	230.0	70.0	60.0	130.0	113.0	NaN	82.0	5.0	...	75.0	NaN	80.0	23.0	7	0.0	1.0	0.0	0.0	0.0
49266	285	2023-03-14 18:00:00+01:00	250.0	70.0	70.0	130.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
60704	313	2023-01-10 08:00:00+01:00	190.0	90.0	100.0	120.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
29283	269	2023-02-20 03:00:00+01:00	230.0	80.0	80.0	110.0	73.0	NaN	63.0	0.0	...	57.0	8.0	93.0	10.0	7	0.0	0.0	0.0	0.0	0.0
9456	240	2023-02-04 00:00:00+01:00	300.0	40.0	40.0	60.0	62.0	NaN	58.0	0.0	...	59.0	8.0	97.0	10.0	7	0.0	0.0	0.0	0.0	0.0
89943	370	2023-02-27 15:00:00+01:00	40.0	60.0	60.0	100.0	57.0	NaN	-38.0	7.0	...	82.0	7.0	50.0	NaN	5	0.0	0.0	0.0	0.0	0.0
58609	312	2023-01-13 01:00:00+01:00	220.0	130.0	120.0	170.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
55103	308	2023-02-15 23:00:00+01:00	180.0	50.0	50.0	80.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
32576	273	2023-01-08 08:00:00+01:00	190.0	50.0	60.0	80.0	75.0	NaN	55.0	5.0	...	70.0	NaN	87.0	NaN	5	0.0	0.0	0.0	0.0	0.0
35926	275	2023-02-26 22:00:00+01:00	30.0	30.0	40.0	50.0	-12.0	NaN	-43.0	0.0	...	75.0	0.0	79.0	NaN	5	0.0	0.0	0.0	0.0	0.0
37712	277	2023-02-11 08:00:00+01:00	280.0	90.0	80.0	110.0	70.0	NaN	58.0	0.0	...	NaN	NaN	92.0	NaN	6	NaN	NaN	NaN	NaN	NaN
55803	308	2023-03-17 03:00:00+01:00	170.0	50.0	40.0	80.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
3002	215	2023-02-05 02:00:00+01:00	270.0	40.0	50.0	70.0	77.0	NaN	57.0	0.0	...	70.0	8.0	87.0	NaN	5	0.0	0.0	0.0	0.0	0.0
62020	313	2023-03-06 04:00:00+01:00	310.0	60.0	60.0	80.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
15147	249	2023-01-02 03:00:00+01:00	240.0	40.0	40.0	60.0	108.0	NaN	105.0	0.0	...	19.0	NaN	98.0	58.0	7	0.0	1.0	0.0	0.0	0.0
65821	316	2023-02-12 13:00:00+01:00	140.0	30.0	20.0	30.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
17129	249	2023-03-25 17:00:00+01:00	270.0	80.0	80.0	120.0	89.0	86.0	61.0	5.0	...	70.0	NaN	82.0	0.0	7	0.0	0.0	0.0	0.0	0.0
58813	312	2023-01-21 13:00:00+01:00	90.0	40.0	40.0	60.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
53228	290	2023-02-27 20:00:00+01:00	30.0	10.0	10.0	10.0	-29.0	NaN	-49.0	0.0	...	70.0	0.0	86.0	NaN	5	0.0	0.0	0.0	0.0	0.0

20 rows × 24 columns

	STN	LON(east)	LAT(north)	ALT(m)	NAME
36	324	4.006	51.596	0.0	Stavenisse
37	330	4.122	51.992	11.9	HoekvanHolland
20	278	6.259	52.435	3.6	Heino
29	311	3.672	51.379	0.0	Hoofdplaat
9	251	5.346	53.392	0.7	HoornTerschelling
46	375	5.707	51.659	22.0	Volkel
2	215	4.437	52.141	-1.1	Voorschoten
7	248	5.174	52.634	0.8	Wijdenes
24	285	6.399	53.575	0.0	Huibertgat
22	280	6.585	53.125	5.2	Eelde
11	258	5.401	52.649	7.3	Houtribdijk
25	286	7.150	53.196	-0.2	NieuwBeerta
15	269	5.520	52.458	-3.7	Lelystad
28	310	3.596	51.442	8.0	Vlissingen
49	391	6.197	51.498	19.5	Arcen
32	315	3.998	51.447	0.0	Hansweert
39	340	4.342	51.449	19.2	Woensdrecht
5	240	4.790	52.318	-3.3	Schiphol
0	209	4.518	52.465	0.0	IJmond
44	356	5.146	51.859	0.7	Herwijnen

STN               0
start_time        0
DD             2499
FH             2500
FF             2499
FX             2500
temp          25920
T10N          87120
TD            25920
SQ            30240
Q             30240
DR            28081
RH            28080
P             47520
VV            45361
N             56160
U             25920
WW            73164
IX                0
M             45415
R             45415
S             45415
O             45415
Y             45415
dtype: int64

STN           0
LON(east)     0
LAT(north)    0
ALT(m)        0
NAME          0
dtype: int64

🏟️ Treinstations

Hier laad ik voorbeelden van de treinstations. En daarna kijk ik of er Null-waarden in de dataset zitten.
Het lijkt erop dat er geen NULL-waarden zijn, De een Null-waarde in code is omdat die Null-waarde is eigenlijk een string die NaN is. en is voor nieuw Amsterdam

display(train_station_locations.sample(sample_size))

display(train_station_locations[train_station_locations['code'].isna()])
display(train_station_locations.isna().sum())

	id	code	uic	name_short	name_medium	name_long	slug	country	type	geo_lat	geo_lng
585	677	ZUE	8503000	Zürich HB	Zürich HB	Zürich HB	zurich-hb	CH	intercitystation	47.378190	8.539200
13	26	AMR	8400050	Alkmaar	Alkmaar	Alkmaar	alkmaar	NL	knooppuntIntercitystation	52.637779	4.739722
561	503	WD	8400702	Woerden	Woerden	Woerden	woerden	NL	knooppuntStoptreinstation	52.084999	4.893611
393	771	NDKP	8821105	N' kempen	Noorderkempen	Noorderkempen	noorderkempen	B	knooppuntIntercitystation	51.358088	4.633977
370	334	MINDEN	8013578	Minden	Minden (Westf)	Minden (Westf)	minden-westf	D	stoptreinstation	52.290138	8.934728
468	792	FSN	8894508	Sint-Nikla	Sint-Niklaas	Sint-Niklaas	sint-niklaas	B	stoptreinstation	51.171944	4.144167
550	524	WZ	8400690	Wezep	Wezep	Wezep	wezep	NL	stoptreinstation	52.454166	6.002500
325	321	LW	8400387	Leeuwarden	Leeuwarden	Leeuwarden	leeuwarden	NL	knooppuntIntercitystation	53.195835	5.792222
93	833	KBB	8015577	Bonn-Beuel	Bonn-Beuel	Bonn-Beuel	bonn-beuel	D	intercitystation	50.738611	7.127778
473	430	SKND	8400574	Sneek N	Sneek Noord	Sneek Noord	sneek-noord	NL	stoptreinstation	53.040932	5.663181
187	171	EM	8010288	Emmerich	Emmerich	Emmerich	emmerich	D	knooppuntStoptreinstation	51.834283	6.255718
214	194	GBR	8400262	Glanerbrug	Glanerbrug	Glanerbrug	glanerbrug	NL	stoptreinstation	52.218613	6.974444
544	517	WT	8400684	Weert	Weert	Weert	weert	NL	intercitystation	51.248611	5.703611
521	849	MEREL	8302593	Venezia	Venezia Santa Lu	Venezia Santa Lucia	venezia-santa-lucia	I	intercitystation	45.441073	12.318855
237	226	HDB	8400293	Hardenberg	Hardenberg	Hardenberg	hardenberg	NL	sneltreinstation	52.572777	6.628611
256	247	HM	8400313	Helmond	Helmond	Helmond	helmond	NL	knooppuntIntercitystation	51.475555	5.661944
61	558	OEYNH	8013585	Bad Oeynh	Bad Oeynhausen	Bad Oeynhausen	bad-oeynhausen	D	intercitystation	52.205461	8.796932
220	846	HOJE	7401318	Göteborg C	Göteborg Central	Göteborg Central	goteborg-central	S	intercitystation	57.708889	11.973333
324	307	LDM	8400386	Leerdam	Leerdam	Leerdam	leerdam	NL	stoptreinstation	51.894722	5.091389
401	634	EOP	8021152	Ochtrup	Ochtrup	Ochtrup	ochtrup	D	stoptreinstation	52.201703	7.182643

id             0
code           1
uic            0
name_short     0
name_medium    0
name_long      0
slug           0
country        0
type           0
geo_lat        0
geo_lng        0
dtype: int64

🔍 Vooronderzoek

Aangezien ik meerdere datasets heb, moet ik eerst even kijken waar de vertragingen zich bevinden.

Hiervoor maak ik een taartdiagram en kijk ik waar de vertragingen zich voordoen.

🚄 Treinverkeer

Bij het treinverkeer is het belangrijk om te weten dat als er 'null' staat, dit nog steeds een betekenis heeft.
Het betekend namelijk dat de event 0 lang is of niet heeft plaatsgevonden.

bins =  list(range(1, 10, 1)) + list(range(10, 70, 10)) + list(range(70, int(Train_trafic_2023['Stop:Departure delay'].max()) +100, 100))

labels = [f'{bins[i]}-{bins[i+1]}' for i in range(len(bins)-1)]

counts = pd.cut(Train_trafic_2023['Stop:Departure delay'], bins=bins, labels=labels, include_lowest=True).dropna().astype(str).value_counts(sort=True)

plt.figure(figsize=(20,6))
plt.barh( counts.index, counts)
plt.title(label='aankomst vertraging in minuten. zonder 0 & nan')
plt.show()

labels = ['0'] + labels

bins = [0] + bins 

counts = pd.cut(Train_trafic_2023['Stop:Departure delay'], bins=bins, labels=labels, include_lowest=True).astype(str).value_counts(sort=True)

plt.figure(figsize=(20,6))
plt.barh( counts.index, counts)
plt.title(label='aankomst vertraging in minuten. zonder 0 & nan')
plt.show()

Uit de grafieken hierboven blijkt dat de meeste treinen op tijd vertrekken van hun locatie. Houd er rekening mee dat zowel '0' als 'nan' aangeven dat er geen vertraging is gemeld.

Als ik de 'nan' en '0' weghaal, blijkt dat de vertraging voornamelijk binnen de 10 minuten ligt.

Maar dit hierboven betreft alleen de vertrekvertraging. Hieronder heb ik de aankomstvertraging op dezelfde manier geanalyseerd.

bins = list(range(1, 10, 1)) + list(range(10, 60, 10)) + list(range(60, int(Train_trafic_2023['Stop:Arrival delay'].max()) +100, 100))

labels = [f'{bins[i]}-{bins[i+1]}' for i in range(len(bins)-1)]

counts = pd.cut(Train_trafic_2023['Stop:Arrival delay'], bins=bins, labels=labels, include_lowest=True).dropna().astype(str).value_counts(sort=True)

plt.figure(figsize=(20,6))
plt.barh( counts.index, counts)
plt.title(label='aankomst vertraging in minuten. zonder 0 & nan')
plt.show()

labels = ['0'] + labels 

bins = [0] + bins

counts = pd.cut(Train_trafic_2023['Stop:Arrival delay'], bins=bins, labels=labels, include_lowest=True).astype(str).value_counts(sort=True)

plt.figure(figsize=(20,6))
plt.barh( counts.index, counts)
plt.title(label='aankomst vertraging in minuten. met 0 & nan')
plt.show()

Als ik naar de vier grafieken kijk, lijken ze goed overeen te komen. Dit betekent dat als een trein meestal te laat aankomt, hij waarschijnlijk ook te laat zal vertrekken.

Met de tabel hieronder ga ik kijken of er op een bepaalde lijn meer vertraging is dan op een andere. Hierbij moet ik wel oppassen en controleren dat ik niet te vaak dezelfde lijn meeneem. Aangezien er is aangegeven dat een RDT-ID meerdere keren kan voorkomen op verschillende dagen, moet ik eerst alles goed bij elkaar zetten. Ik neem het gemiddelde omdat ik de gemiddelde vertraging van de ritten op een dag wil zien.

grouped = Train_trafic_2023.groupby(['Service:Date', 'Service:RDT-ID'])['Stop:Arrival delay'].mean()

df_avarage_delay_day = grouped.reset_index(name='Stop:Arrival delay')

grouped = df_avarage_delay_day.groupby(['Service:Date'])['Stop:Arrival delay'].mean()

df_avarage_delay_day = grouped.reset_index(name='Stop:Arrival delay')

display(df_avarage_delay_day)

plt.figure(figsize=(20,6))
plt.bar( df_avarage_delay_day['Service:Date'], df_avarage_delay_day['Stop:Arrival delay'])
plt.title(label='avarage aankomst vertraging door de maanden.')
plt.show()

	Service:Date	Stop:Arrival delay
0	2023-01-01	0.680368
1	2023-01-02	0.596654
2	2023-01-03	0.565669
3	2023-01-04	0.781114
4	2023-01-05	0.673781
...	...	...
85	2023-03-27	0.799672
86	2023-03-28	0.964126
87	2023-03-29	0.642580
88	2023-03-30	0.710096
89	2023-03-31	0.702228

90 rows × 2 columns

In de grafiek hierboven is te zien dat de gemiddelde aankomstvertraging varieert gedurende de maanden.
Dit kan belangrijk zijn voor het testen van de invloed van weersomstandigheden.

🌥️ Het weer

Er zijn veel kolommen in de weergegevens en ik ga eerst bepalen wat nodig is en wat niet. Daarvoor ga ik het gemiddelde berekenen en kijken hoe dicht dat bij het antwoord ligt. Als het dichtbij genoeg is, kan ik concluderen dat het gemiddelde representatief is voor heel Nederland. Als dat niet het geval is, moet ik even nadenken over wat ik dan ga doen.

grouped = weather_2023.groupby(by=['start_time']).mean()

db_weather_single_day = grouped.reset_index()

db_weather_single_day = db_weather_single_day.drop(columns=['STN'])

display(db_weather_single_day.head(24))

df_first_date = weather_2023[weather_2023['start_time'] == weather_2023['start_time'].min()]

display(df_first_date[df_first_date['temp'].notna()].sort_values('temp'))

	start_time	DD	FH	FF	FX	temp	T10N	TD	SQ	Q	...	VV	N	U	WW	IX	M	R	S	O	Y
0	2023-01-01 00:00:00+01:00	215.555556	110.666667	110.000000	172.222222	140.617647	NaN	80.382353	0.00000	0.00000	...	73.28	8.00	67.558824	36.062500	6.152174	0.0	0.56	0.0	0.0	0.0
1	2023-01-01 01:00:00+01:00	215.333333	108.000000	106.444444	167.333333	138.029412	NaN	80.764706	0.00000	0.00000	...	74.28	8.00	68.941176	49.187500	6.152174	0.0	0.64	0.0	0.0	0.0
2	2023-01-01 02:00:00+01:00	220.222222	107.333333	106.888889	168.000000	135.529412	NaN	79.941176	0.00000	0.00000	...	74.12	8.00	69.764706	30.631579	6.282609	0.0	0.48	0.0	0.0	0.0
3	2023-01-01 03:00:00+01:00	223.111111	104.000000	103.333333	165.111111	131.323529	NaN	80.911765	0.00000	0.00000	...	73.20	7.95	72.147059	34.062500	6.152174	0.0	0.64	0.0	0.0	0.0
4	2023-01-01 04:00:00+01:00	224.444444	99.111111	99.111111	158.444444	125.970588	NaN	79.970588	0.00000	0.00000	...	73.88	8.00	73.911765	20.833333	5.978261	0.0	0.40	0.0	0.0	0.0
5	2023-01-01 05:00:00+01:00	224.666667	95.333333	94.222222	148.444444	122.558824	114.294118	76.294118	0.00000	0.00000	...	73.08	8.00	73.529412	11.791667	6.500000	0.0	0.20	0.0	0.0	0.0
6	2023-01-01 06:00:00+01:00	224.666667	91.333333	89.555556	140.222222	117.941176	NaN	75.352941	0.00000	0.00000	...	70.80	8.00	75.117647	45.666667	5.978261	0.0	0.40	0.0	0.0	0.0
7	2023-01-01 07:00:00+01:00	222.666667	87.111111	83.777778	134.444444	113.647059	NaN	75.970588	0.00000	0.25000	...	69.88	8.00	77.529412	36.153846	6.021739	0.0	0.40	0.0	0.0	0.0
8	2023-01-01 08:00:00+01:00	219.111111	81.777778	81.555556	127.555556	112.441176	NaN	78.558824	0.18750	5.62500	...	69.40	8.00	79.558824	22.785714	6.065217	0.0	0.24	0.0	0.0	0.0
9	2023-01-01 09:00:00+01:00	215.111111	76.666667	76.444444	119.555556	114.911765	NaN	81.676471	2.62500	25.12500	...	69.56	8.00	79.852941	28.166667	5.717391	0.0	0.16	0.0	0.0	0.0
10	2023-01-01 10:00:00+01:00	216.222222	78.444444	80.222222	120.000000	118.764706	NaN	83.588235	3.25000	40.75000	...	70.64	7.90	78.823529	16.500000	5.630435	0.0	0.08	0.0	0.0	0.0
11	2023-01-01 11:00:00+01:00	218.444444	78.666667	77.777778	122.888889	120.529412	102.852941	84.882353	1.46875	42.56250	...	71.48	7.95	78.441176	10.666667	5.978261	0.0	0.08	0.0	0.0	0.0
12	2023-01-01 12:00:00+01:00	218.666667	77.555556	74.222222	119.555556	121.382353	NaN	84.088235	2.00000	38.59375	...	71.24	8.00	77.735294	14.333333	5.586957	0.0	0.04	0.0	0.0	0.0
13	2023-01-01 13:00:00+01:00	212.222222	68.000000	63.777778	107.555556	119.647059	NaN	83.529412	1.06250	23.28125	...	71.04	8.00	78.147059	10.000000	5.543478	0.0	0.00	0.0	0.0	0.0
14	2023-01-01 14:00:00+01:00	201.555556	59.111111	54.888889	92.666667	116.352941	NaN	83.352941	1.15625	10.03125	...	70.52	8.00	79.852941	25.750000	5.630435	0.0	0.04	0.0	0.0	0.0
15	2023-01-01 15:00:00+01:00	190.000000	52.000000	49.777778	82.000000	113.617647	NaN	83.647059	0.00000	1.18750	...	69.88	8.00	81.441176	56.285714	5.760870	0.0	0.24	0.0	0.0	0.0
16	2023-01-01 16:00:00+01:00	182.000000	48.888889	49.777778	76.666667	111.058824	NaN	86.264706	0.00000	0.00000	...	68.64	8.00	84.294118	64.555556	5.847826	0.0	0.36	0.0	0.0	0.0
17	2023-01-01 17:00:00+01:00	176.000000	49.111111	49.555556	78.222222	108.735294	100.294118	89.058824	0.00000	0.00000	...	67.32	8.00	87.323529	61.000000	6.108696	0.0	0.60	0.0	0.0	0.0
18	2023-01-01 18:00:00+01:00	172.000000	49.333333	49.555556	78.444444	105.764706	NaN	92.147059	0.00000	0.00000	...	65.08	8.00	90.970588	53.947368	6.282609	0.0	0.76	0.0	0.0	0.0
19	2023-01-01 19:00:00+01:00	173.555556	48.222222	48.222222	78.222222	104.058824	NaN	93.617647	0.00000	0.00000	...	62.68	8.00	92.823529	62.083333	6.500000	0.0	0.96	0.0	0.0	0.0
20	2023-01-01 20:00:00+01:00	187.333333	41.333333	40.666667	76.222222	101.941176	NaN	93.500000	0.00000	0.00000	...	60.88	8.00	94.117647	59.280000	6.543478	0.0	1.00	0.0	0.0	0.0
21	2023-01-01 21:00:00+01:00	167.555556	40.444444	40.000000	73.777778	102.264706	NaN	94.764706	0.00000	0.00000	...	61.36	8.00	94.705882	52.920000	6.543478	0.0	1.00	0.0	0.0	0.0
22	2023-01-01 22:00:00+01:00	172.000000	41.777778	43.555556	70.000000	105.029412	NaN	95.470588	0.00000	0.00000	...	60.68	8.00	93.529412	46.590909	6.413043	0.0	0.88	0.0	0.0	0.0
23	2023-01-01 23:00:00+01:00	171.777778	43.777778	42.444444	72.444444	107.970588	95.647059	96.911765	0.00000	0.00000	...	59.48	8.00	92.647059	42.840000	6.543478	0.0	0.84	0.0	0.0	0.0

24 rows × 23 columns

	STN	start_time	DD	FH	FF	FX	temp	T10N	TD	SQ	...	VV	N	U	WW	IX	M	R	S	O	Y
17280	251	2023-01-01 00:00:00+01:00	220.0	140.0	130.0	190.0	108.0	NaN	88.0	0.0	...	58.0	NaN	87.0	10.0	7	0.0	0.0	0.0	0.0	0.0
6480	235	2023-01-01 00:00:00+01:00	220.0	120.0	120.0	180.0	112.0	NaN	95.0	0.0	...	56.0	8.0	89.0	52.0	7	0.0	1.0	0.0	0.0	0.0
10800	242	2023-01-01 00:00:00+01:00	220.0	160.0	150.0	210.0	114.0	NaN	91.0	NaN	...	57.0	8.0	86.0	10.0	7	0.0	0.0	0.0	0.0	0.0
25920	267	2023-01-01 00:00:00+01:00	220.0	130.0	130.0	200.0	121.0	NaN	89.0	0.0	...	64.0	8.0	81.0	57.0	7	0.0	1.0	0.0	0.0	0.0
19440	257	2023-01-01 00:00:00+01:00	NaN	NaN	NaN	NaN	123.0	NaN	83.0	0.0	...	NaN	NaN	76.0	NaN	6	NaN	NaN	NaN	NaN	NaN
49680	286	2023-01-01 00:00:00+01:00	210.0	110.0	100.0	170.0	128.0	NaN	98.0	0.0	...	NaN	NaN	82.0	NaN	6	NaN	NaN	NaN	NaN	NaN
36720	277	2023-01-01 00:00:00+01:00	220.0	120.0	130.0	200.0	130.0	NaN	92.0	0.0	...	NaN	NaN	78.0	NaN	6	NaN	NaN	NaN	NaN	NaN
32400	273	2023-01-01 00:00:00+01:00	220.0	80.0	70.0	130.0	131.0	NaN	93.0	0.0	...	64.0	NaN	77.0	23.0	7	0.0	1.0	0.0	0.0	0.0
56160	310	2023-01-01 00:00:00+01:00	200.0	140.0	140.0	180.0	131.0	NaN	79.0	0.0	...	75.0	8.0	70.0	NaN	5	0.0	0.0	0.0	0.0	0.0
30240	270	2023-01-01 00:00:00+01:00	230.0	140.0	140.0	230.0	132.0	NaN	94.0	0.0	...	74.0	8.0	78.0	23.0	7	0.0	1.0	0.0	0.0	0.0
43200	280	2023-01-01 00:00:00+01:00	210.0	110.0	100.0	190.0	133.0	NaN	94.0	0.0	...	71.0	8.0	77.0	23.0	7	0.0	1.0	0.0	0.0	0.0
73440	330	2023-01-01 00:00:00+01:00	210.0	130.0	120.0	210.0	135.0	NaN	78.0	0.0	...	NaN	NaN	68.0	NaN	6	NaN	NaN	NaN	NaN	NaN
15120	249	2023-01-01 00:00:00+01:00	230.0	140.0	130.0	230.0	136.0	NaN	98.0	0.0	...	63.0	NaN	77.0	81.0	7	0.0	1.0	0.0	0.0	0.0
8640	240	2023-01-01 00:00:00+01:00	220.0	110.0	110.0	150.0	138.0	NaN	92.0	0.0	...	69.0	8.0	73.0	81.0	7	0.0	1.0	0.0	0.0	0.0
41040	279	2023-01-01 00:00:00+01:00	220.0	90.0	90.0	150.0	139.0	NaN	92.0	0.0	...	75.0	8.0	73.0	51.0	7	0.0	1.0	0.0	0.0	0.0
28080	269	2023-01-01 00:00:00+01:00	220.0	100.0	110.0	170.0	141.0	NaN	90.0	0.0	...	67.0	8.0	71.0	23.0	7	0.0	1.0	0.0	0.0	0.0
2160	215	2023-01-01 00:00:00+01:00	220.0	110.0	130.0	200.0	142.0	NaN	84.0	0.0	...	66.0	8.0	68.0	NaN	5	0.0	0.0	0.0	0.0	0.0
79920	344	2023-01-01 00:00:00+01:00	210.0	100.0	110.0	180.0	143.0	NaN	88.0	0.0	...	79.0	8.0	69.0	NaN	5	0.0	0.0	0.0	0.0	0.0
66960	319	2023-01-01 00:00:00+01:00	210.0	100.0	90.0	180.0	145.0	NaN	74.0	0.0	...	81.0	NaN	62.0	NaN	5	0.0	0.0	0.0	0.0	0.0
69120	323	2023-01-01 00:00:00+01:00	210.0	90.0	90.0	160.0	145.0	NaN	80.0	0.0	...	NaN	NaN	65.0	NaN	6	NaN	NaN	NaN	NaN	NaN
95040	380	2023-01-01 00:00:00+01:00	210.0	100.0	90.0	150.0	147.0	NaN	42.0	0.0	...	83.0	8.0	49.0	NaN	5	0.0	0.0	0.0	0.0	0.0
82080	348	2023-01-01 00:00:00+01:00	220.0	120.0	130.0	180.0	147.0	NaN	82.0	0.0	...	76.0	NaN	64.0	NaN	5	0.0	0.0	0.0	0.0	0.0
75600	340	2023-01-01 00:00:00+01:00	220.0	70.0	80.0	140.0	148.0	NaN	72.0	NaN	...	83.0	8.0	60.0	23.0	7	0.0	1.0	0.0	0.0	0.0
38880	278	2023-01-01 00:00:00+01:00	220.0	70.0	60.0	130.0	150.0	NaN	83.0	0.0	...	NaN	NaN	64.0	NaN	6	NaN	NaN	NaN	NaN	NaN
86400	356	2023-01-01 00:00:00+01:00	220.0	110.0	110.0	170.0	151.0	NaN	80.0	0.0	...	NaN	NaN	62.0	NaN	6	NaN	NaN	NaN	NaN	NaN
51840	290	2023-01-01 00:00:00+01:00	210.0	80.0	70.0	150.0	152.0	NaN	72.0	0.0	...	78.0	8.0	58.0	23.0	7	0.0	1.0	0.0	0.0	0.0
34560	275	2023-01-01 00:00:00+01:00	220.0	100.0	90.0	150.0	153.0	NaN	76.0	0.0	...	82.0	8.0	59.0	23.0	7	0.0	1.0	0.0	0.0	0.0
45360	283	2023-01-01 00:00:00+01:00	210.0	70.0	70.0	130.0	154.0	NaN	72.0	0.0	...	NaN	NaN	58.0	NaN	6	NaN	NaN	NaN	NaN	NaN
84240	350	2023-01-01 00:00:00+01:00	210.0	90.0	90.0	150.0	154.0	NaN	73.0	0.0	...	83.0	8.0	58.0	51.0	7	0.0	1.0	0.0	0.0	0.0
23760	260	2023-01-01 00:00:00+01:00	220.0	80.0	100.0	160.0	155.0	NaN	75.0	0.0	...	81.0	8.0	58.0	NaN	5	0.0	0.0	0.0	0.0	0.0
92880	377	2023-01-01 00:00:00+01:00	210.0	90.0	80.0	160.0	155.0	NaN	53.0	0.0	...	84.0	8.0	50.0	NaN	5	0.0	0.0	0.0	0.0	0.0
90720	375	2023-01-01 00:00:00+01:00	220.0	110.0	110.0	160.0	162.0	NaN	66.0	0.0	...	81.0	8.0	52.0	NaN	5	0.0	0.0	0.0	0.0	0.0
88560	370	2023-01-01 00:00:00+01:00	220.0	110.0	120.0	190.0	163.0	NaN	59.0	0.0	...	82.0	8.0	49.0	23.0	7	0.0	1.0	0.0	0.0	0.0
97200	391	2023-01-01 00:00:00+01:00	220.0	50.0	50.0	120.0	163.0	NaN	56.0	0.0	...	NaN	NaN	49.0	NaN	6	NaN	NaN	NaN	NaN	NaN

34 rows × 24 columns

Het lijkt erop dat er met de temperatuur een verschil van ongeveer 4 graden is tussen de hoogste en laagste temperatuur. Hoewel een verschil van 4 graden niet veel lijkt, kan het ervoor zorgen dat mijn data niet meer klopt. Het lijkt er dus op dat ik toch het verschil tussen de stations moet vinden en het dichtstbijzijnde weerstation moet bepalen. Daarbij moet ik ook rekening houden met het tijdstip van de meting.

🏟️ Treinstations

Allereerst moet ik bepalen welke treinstations ik nodig heb voor mijn datasets. Dat is wat ik hieronder doe. Ik neem de stationscodes van het treinverkeer en de codes van de treinstations in de locatiedataset en kijk welke er overeenkomen.

train_station_locations = train_station_locations[train_station_locations['code'].isin(Train_trafic_2023['Stop:Station code'])]

display(train_station_locations)

	id	code	uic	name_short	name_medium	name_long	slug	country	type	geo_lat	geo_lng
0	266	HT	8400319	Den Bosch	's-Hertogenbosch	's-Hertogenbosch	s-hertogenbosch	NL	knooppuntIntercitystation	51.690480	5.293620
1	269	HTO	8400320	Dn Bosch O	's-Hertogenb. O.	's-Hertogenbosch Oost	s-hertogenbosch-oost	NL	stoptreinstation	51.700554	5.318333
2	227	HDE	8400388	't Harde	't Harde	't Harde	t-harde	NL	stoptreinstation	52.409168	5.893611
3	8	AHBF	8015345	Aachen	Aachen Hbf	Aachen Hbf	aachen-hbf	D	knooppuntIntercitystation	50.767800	6.091499
4	818	AW	8015199	Aachen W	Aachen West	Aachen West	aachen-west	D	stoptreinstation	50.780360	6.070715
...	...	...	...	...	...	...	...	...	...	...	...
585	677	ZUE	8503000	Zürich HB	Zürich HB	Zürich HB	zurich-hb	CH	intercitystation	47.378190	8.539200
586	534	ZP	8400744	Zutphen	Zutphen	Zutphen	zutphen	NL	knooppuntIntercitystation	52.145279	6.194167
587	541	ZWD	8400752	Zwijndrcht	Zwijndrecht	Zwijndrecht	zwijndrecht	NL	stoptreinstation	51.814999	4.641667
589	532	ZL	8400747	Zwolle	Zwolle	Zwolle	zwolle	NL	knooppuntIntercitystation	52.504723	6.091944
590	815	ZLSH	8400746	Stadshagen	Stadshagen	Zwolle Stadshagen	zwolle-stadshagen	NL	stoptreinstation	52.527640	6.051460

546 rows × 11 columns

Met deze informatie kan ik mijn treinstations coördinaten geven, zodat ik het dichtstbijzijnde weerstation kan gebruiken.

Maar om dat te testen, voer ik hier eerst een test uit gebaseerd op één meetmoment.

Hier gebruik ik de BallTree-functie van de 'neighbors'-functies van sklearn. Hiermee probeer ik de temperaturen van de stations te ontdekken.

De 'haversine'-metric is een goede metric omdat het de afstand tussen twee locaties kan berekenen.

Ik moet alleen voor de 'haversine'-metric de invoerinformatie omzetten naar de radiaalvorm, gelukkig heeft numpy daar een goede functie voor.

Uiteindelijk doe ik dit om een zo dicht mogelijk weerstation te krijgen. Dit heb ik nodig om te bepalen of de vertragingen in het treinverkeer een verband hebben met het weer.

from sklearn.neighbors import BallTree

df_mergid_first_date_stations = df_first_date.merge(weather_station_location, on='STN', how='inner')

df_mergid_first_date_stations[['LAT(north)_radians', 'LON(east)_radians']] = np.radians(df_mergid_first_date_stations[['LAT(north)', 'LON(east)']])
train_station_locations[['geo_lat_radians', 'geo_lng_radians']] = np.radians(train_station_locations[['geo_lat', 'geo_lng']])

tree = BallTree(df_mergid_first_date_stations[['LAT(north)_radians', 'LON(east)_radians']].dropna().values, metric='haversine')

indices = tree.query(train_station_locations[['geo_lat_radians', 'geo_lng_radians']].dropna().values, return_distance=False, k=1)

nearest_neighbors_weather_station_location = df_mergid_first_date_stations.iloc[indices.flatten()]

train_station_locations.reset_index(drop=True, inplace=True)
nearest_neighbors_weather_station_location.reset_index(drop=True, inplace=True)

combined_df = pd.concat([nearest_neighbors_weather_station_location, train_station_locations], axis=1)

display(combined_df.drop(columns=['geo_lat_radians', 'geo_lng_radians', 'LAT(north)_radians', 'LON(east)_radians', 'name_short', 'name_medium', 'name_long', 'id', 'code', 'uic']))

	STN	start_time	DD	FH	FF	FX	temp	T10N	TD	SQ	...	Y	LON(east)	LAT(north)	ALT(m)	NAME	slug	country	type	geo_lat	geo_lng
0	356	2023-01-01 00:00:00+01:00	220.0	110.0	110.0	170.0	151.0	NaN	80.0	0.0	...	NaN	5.146	51.859	0.7	Herwijnen	s-hertogenbosch	NL	knooppuntIntercitystation	51.690480	5.293620
1	356	2023-01-01 00:00:00+01:00	220.0	110.0	110.0	170.0	151.0	NaN	80.0	0.0	...	NaN	5.146	51.859	0.7	Herwijnen	s-hertogenbosch-oost	NL	stoptreinstation	51.700554	5.318333
2	278	2023-01-01 00:00:00+01:00	220.0	70.0	60.0	130.0	150.0	NaN	83.0	0.0	...	NaN	6.259	52.435	3.6	Heino	t-harde	NL	stoptreinstation	52.409168	5.893611
3	380	2023-01-01 00:00:00+01:00	210.0	100.0	90.0	150.0	147.0	NaN	42.0	0.0	...	0.0	5.762	50.906	114.3	Maastricht	aachen-hbf	D	knooppuntIntercitystation	50.767800	6.091499
4	380	2023-01-01 00:00:00+01:00	210.0	100.0	90.0	150.0	147.0	NaN	42.0	0.0	...	0.0	5.762	50.906	114.3	Maastricht	aachen-west	D	stoptreinstation	50.780360	6.070715
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
541	380	2023-01-01 00:00:00+01:00	210.0	100.0	90.0	150.0	147.0	NaN	42.0	0.0	...	0.0	5.762	50.906	114.3	Maastricht	zurich-hb	CH	intercitystation	47.378190	8.539200
542	275	2023-01-01 00:00:00+01:00	220.0	100.0	90.0	150.0	153.0	NaN	76.0	0.0	...	0.0	5.873	52.056	48.2	Deelen	zutphen	NL	knooppuntIntercitystation	52.145279	6.194167
543	344	2023-01-01 00:00:00+01:00	210.0	100.0	110.0	180.0	143.0	NaN	88.0	0.0	...	0.0	4.447	51.962	-4.3	Rotterdam	zwijndrecht	NL	stoptreinstation	51.814999	4.641667
544	278	2023-01-01 00:00:00+01:00	220.0	70.0	60.0	130.0	150.0	NaN	83.0	0.0	...	NaN	6.259	52.435	3.6	Heino	zwolle	NL	knooppuntIntercitystation	52.504723	6.091944
545	278	2023-01-01 00:00:00+01:00	220.0	70.0	60.0	130.0	150.0	NaN	83.0	0.0	...	NaN	6.259	52.435	3.6	Heino	zwolle-stadshagen	NL	stoptreinstation	52.527640	6.051460

546 rows × 33 columns

Uit de test hierboven kan ik concluderen dat ik met de weerstations en treinstationlocaties redelijk goed het dichtstbijzijnde weerstation van een treinstation kan vinden.

🎛️ Correlatie bekijken

Hier ga ik de correlatieplots van de datasets weergeven.

🚄 Treinverkeer

Dit is de correlatieplot van het treinverkeer:

corr = Train_trafic_2023.select_dtypes(include=['int32', 'int64', 'float64']).corr()

mask = np.tril(corr)

sns.heatmap(corr, mask=mask, cbar=True, annot=True, fmt=".1f")

<Axes: >

Het lijkt erop dat als ik de correlatieplot van de afzonderlijke treininformatie maak, er een lage correlatie tussen de gegevens is. Hieronder ga ik een pairplot maken om te zien of ik er meer uit kan halen.

plt.figure(figsize=(20,20))
sns.pairplot(Train_trafic_2023.sample(frac=0.001), hue='Service:Company')

<seaborn.axisgrid.PairGrid at 0x167f56ee510>

<Figure size 2000x2000 with 0 Axes>

Hoewel ik de hoeveelheid data heb moeten verminderen, heb ik toch een goed gevoel over de inhoud van mijn dataset.

🌥️ het weer

plt.subplots(figsize=(15, 15))

corr = weather_2023.select_dtypes(include=['int32', 'int64', 'float64']).corr()
mask = np.tril(corr);
sns.heatmap(corr, mask=mask, cbar=True, annot=True, fmt=".1f",)

<Axes: >

Het lijkt erop dat verschillende weerelementen een sterke relatie met elkaar hebben. Het gaat niet alleen om de hoeveelheid, maar ook om de elementen die naar verwachting met elkaar in verband staan. Bijvoorbeeld, de temperatuur en de minimale temperatuur elke 6 uur vertonen een sterke correlatie.

plt.figure(figsize=(20,20))
display(weather_2023.sample(frac=0.001))
sns.pairplot(weather_2023.sample(frac=0.001), hue='STN')

	STN	start_time	DD	FH	FF	FX	temp	T10N	TD	SQ	...	VV	N	U	WW	IX	M	R	S	O	Y
88492	356	2023-03-29 04:00:00+02:00	190.0	30.0	30.0	50.0	72.0	NaN	34.0	0.0	...	NaN	NaN	76.0	NaN	6	NaN	NaN	NaN	NaN	NaN
87135	356	2023-01-31 15:00:00+01:00	270.0	70.0	70.0	110.0	83.0	NaN	58.0	3.0	...	NaN	NaN	84.0	NaN	6	NaN	NaN	NaN	NaN	NaN
3988	215	2023-03-18 04:00:00+01:00	200.0	30.0	20.0	60.0	93.0	NaN	84.0	0.0	...	60.0	8.0	94.0	23.0	7	0.0	1.0	0.0	0.0	0.0
93522	377	2023-01-27 18:00:00+01:00	20.0	30.0	30.0	50.0	40.0	NaN	5.0	0.0	...	75.0	8.0	77.0	NaN	5	0.0	0.0	0.0	0.0	0.0
1209	209	2023-02-20 09:00:00+01:00	240.0	120.0	120.0	160.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
62738	315	2023-01-05 02:00:00+01:00	260.0	90.0	90.0	130.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
54643	308	2023-01-27 19:00:00+01:00	NaN	NaN	NaN	NaN	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
86680	356	2023-01-12 16:00:00+01:00	230.0	90.0	100.0	150.0	107.0	NaN	99.0	0.0	...	NaN	NaN	94.0	NaN	6	NaN	NaN	NaN	NaN	NaN
27661	267	2023-03-14 13:00:00+01:00	310.0	100.0	100.0	130.0	58.0	NaN	-3.0	9.0	...	76.0	7.0	64.0	NaN	5	0.0	0.0	0.0	0.0	0.0
1600	209	2023-03-08 16:00:00+01:00	90.0	70.0	60.0	90.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN

99 rows × 24 columns

<seaborn.axisgrid.PairGrid at 0x1682db41d50>

<Figure size 2000x2000 with 0 Axes>

🏟️ trein station

corr = train_station_locations.select_dtypes(include=['int32', 'int64', 'float64']).corr()
mask = np.tril(corr)
sns.heatmap(corr, mask=mask, cbar=True, annot=True, fmt=".1f")

<Axes: >

Dit was meer om te kijken of ik iets terug zou krijgen, maar het geeft een goed voorbeeld van de gegevens.

🛠️ Preprocessing

Ik ga onderzoeken of ik de datasets kan samenvoegen, zodat ik één volledige dataset krijg met mijn informatie.
Uiteraard moet de data op de juiste plek staan.

🚄&🏟️ Treinverkeer samenvoegen met treinstation

Hieronder ga ik al het treinverkeer voorzien van hun treinstationgegevens, zodat ik het later kan samenvoegen met het juiste weerstation.
Hiervoor ga ik een nieuwe variabele aanmaken, zodat ik het meerdere keren kan uitvoeren zonder dat het de data hierboven verstoort.
De samenvoeging zal worden gedaan door de twee datasets te joinen met de merge-functie van Panda's.
De join zal alle gegevens van beide tabellen bevatten. Ik gebruik in de merge de stationscode die in beide tabellen voorkomt. Hoewel ze anders worden genoemd, zijn het allebei de stationscodes.
Dit zou niet het moeilijkste moeten zijn, aangezien de datasets van dezelfde provider komen.

# train_station_locations , Train_trafic_2023
display(
   train_station_locations[train_station_locations['name_short']== "Nw A'dam"] 
)
train_trafic_with_locations = Train_trafic_2023.merge(right=train_station_locations, how='outer', left_on='Stop:Station code', right_on='code')
display(train_trafic_with_locations)

	id	code	uic	name_short	name_medium	name_long	slug	country	type	geo_lat	geo_lng	geo_lat_radians	geo_lng_radians
351	349	NaN	8400454	Nw A'dam	Nw. Amsterdam	Nieuw Amsterdam	nieuw-amsterdam	NL	sneltreinstation	52.718613	6.848611	0.920113	0.119531

	Service:RDT-ID	Service:Date	Service:Type	Service:Company	Service:Train number	Service:Completely cancelled	Service:Partly cancelled	Service:Maximum delay	Stop:RDT-ID	Stop:Station code	...	name_short	name_medium	name_long	slug	country	type	geo_lat	geo_lng	geo_lat_radians	geo_lng_radians
0	10368001	2023-01-01	sprinter	ns	4026	False	False	0	93070750	AC	...	Abcoude	Abcoude	Abcoude	abcoude	NL	stoptreinstation	52.278500	4.977000	0.912432	0.086865
1	10368079	2023-01-01	sprinter	ns	4037	False	False	0	93071464	AC	...	Abcoude	Abcoude	Abcoude	abcoude	NL	stoptreinstation	52.278500	4.977000	0.912432	0.086865
2	10368089	2023-01-01	sprinter	ns	4035	False	False	0	93071556	AC	...	Abcoude	Abcoude	Abcoude	abcoude	NL	stoptreinstation	52.278500	4.977000	0.912432	0.086865
3	10368125	2023-01-01	sprinter	ns	4030	False	False	0	93071943	AC	...	Abcoude	Abcoude	Abcoude	abcoude	NL	stoptreinstation	52.278500	4.977000	0.912432	0.086865
4	10368142	2023-01-01	sprinter	ns	4028	False	False	0	93072094	AC	...	Abcoude	Abcoude	Abcoude	abcoude	NL	stoptreinstation	52.278500	4.977000	0.912432	0.086865
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
5142245	10915642	2023-03-31	stoptrein	blauwnet	8092	False	False	0	98208124	NaN	...	Nw A'dam	Nw. Amsterdam	Nieuw Amsterdam	nieuw-amsterdam	NL	sneltreinstation	52.718613	6.848611	0.920113	0.119531
5142246	10915677	2023-03-31	sneltrein	blauwnet	3879	False	False	0	98208425	NaN	...	Nw A'dam	Nw. Amsterdam	Nieuw Amsterdam	nieuw-amsterdam	NL	sneltreinstation	52.718613	6.848611	0.920113	0.119531
5142247	10915796	2023-03-31	stoptrein	blauwnet	8081	False	False	0	98209348	NaN	...	Nw A'dam	Nw. Amsterdam	Nieuw Amsterdam	nieuw-amsterdam	NL	sneltreinstation	52.718613	6.848611	0.920113	0.119531
5142248	10915874	2023-03-31	stoptrein	blauwnet	8096	False	False	0	98209890	NaN	...	Nw A'dam	Nw. Amsterdam	Nieuw Amsterdam	nieuw-amsterdam	NL	sneltreinstation	52.718613	6.848611	0.920113	0.119531
5142249	10915882	2023-03-31	stoptrein	blauwnet	8083	False	False	0	98209958	NaN	...	Nw A'dam	Nw. Amsterdam	Nieuw Amsterdam	nieuw-amsterdam	NL	sneltreinstation	52.718613	6.848611	0.920113	0.119531

5142250 rows × 30 columns

Voor deze verbinding heb ik ontdekt dat het station 'New Amsterdam' een stationscode heeft die 'NaN' is, en niet de variabele NaN.
Maar na deze ontdekking verliep de verbinding vlekkeloos.

🌥️&🚄🏟️ Het weer samenvoegen met het al samengevoegde treinverkeer en stations

Voor deze taak zal ik waarschijnlijk machine learning moeten gebruiken om het juiste weerstation aan een treinstation te koppelen.
Dit komt omdat een weerstation ook in de loop van de tijd kan veranderen, wat invloed heeft op het weer.

Als eerste stap ga ik kijken welk weerstation het dichtst bij welk treinstation ligt.

from sklearn.neighbors import BallTree
import numpy as np

# weather_station_location, train_trafic_with_locations, weather_2023

# Weather_With_location_time = weather_2023.merge(weather_station_location, on='STN', how='inner')

tree = BallTree(np.radians(weather_station_location[['LAT(north)', 'LON(east)']]), metric='haversine')
distances, indices = tree.query(np.radians(train_trafic_with_locations[['geo_lat', 'geo_lng']].dropna()), dualtree=True)

display(distances)
display(indices)
display(weather_2023)

array([[0.00211167],
       [0.00211167],
       [0.00211167],
       ...,
       [0.00295339],
       [0.00295339],
       [0.00295339]])

array([[ 5],
       [ 5],
       [ 5],
       ...,
       [21],
       [21],
       [21]], dtype=int64)

	STN	start_time	DD	FH	FF	FX	temp	T10N	TD	SQ	...	VV	N	U	WW	IX	M	R	S	O	Y
0	209	2023-01-01 00:00:00+01:00	220.0	140.0	150.0	170.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
1	209	2023-01-01 01:00:00+01:00	220.0	130.0	130.0	170.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
2	209	2023-01-01 02:00:00+01:00	240.0	140.0	130.0	200.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
3	209	2023-01-01 03:00:00+01:00	220.0	120.0	130.0	170.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
4	209	2023-01-01 04:00:00+01:00	230.0	130.0	140.0	170.0	NaN	NaN	NaN	NaN	...	NaN	NaN	NaN	NaN	6	NaN	NaN	NaN	NaN	NaN
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
99355	391	2023-03-31 19:00:00+02:00	210.0	50.0	40.0	120.0	99.0	NaN	88.0	0.0	...	NaN	NaN	93.0	NaN	6	NaN	NaN	NaN	NaN	NaN
99356	391	2023-03-31 20:00:00+02:00	220.0	40.0	50.0	100.0	99.0	NaN	86.0	0.0	...	NaN	NaN	91.0	NaN	6	NaN	NaN	NaN	NaN	NaN
99357	391	2023-03-31 21:00:00+02:00	210.0	40.0	40.0	110.0	97.0	NaN	87.0	0.0	...	NaN	NaN	93.0	NaN	6	NaN	NaN	NaN	NaN	NaN
99358	391	2023-03-31 22:00:00+02:00	220.0	40.0	40.0	90.0	95.0	NaN	88.0	0.0	...	NaN	NaN	95.0	NaN	6	NaN	NaN	NaN	NaN	NaN
99359	391	2023-03-31 23:00:00+02:00	220.0	50.0	50.0	110.0	94.0	96.0	88.0	0.0	...	NaN	NaN	95.0	NaN	6	NaN	NaN	NaN	NaN	NaN

99360 rows × 24 columns

Hierboven heb ik ervoor gezorgd dat ik weet welk weerstation het dichtst bij een treinstation staat.
Deze informatie heb ik hieronder nodig om een dataset te creëren die alle gegevens van het weerstation bij de treingegevens voegt.

Hieronder ga ik de data samenvoegen zodat ik uiteindelijk één dataset heb waar alles in staat dat ik nodig zou moeten hebben.
Dit doe ik door de dichtstbijzijnde tijd bij zijn dichtstbijzijnde weerstation te zoeken. Ik heb ervoor gekozen om voor het testen minder data te gebruiken, zodat ik sneller feedback krijg van mijn outputs.

En omdat het een grote dataset is ga ik hem opslaan naar een csv om hem dan te laten kiezen tussen de dataframe maken of de dataframe inladen van geschiedenis.
uitijndelijk do ik dit om het te versnellen en ruimte in mijn RAM geheugen op te schonen.

CSV_save_location_full_data_weather_and_train = datafolder + 'Processed data/Weather and train station data.csv'

if(calculate_large_data): 
    combined_data = []

    for i, train_station in train_trafic_with_locations.iterrows():
        if(indices.size > i):
            nearest_weather_station = weather_station_location.iloc[indices[i][0]]
            nearest_weather_station = weather_2023[weather_2023['STN'] == nearest_weather_station.STN]
            if(nearest_weather_station['start_time'].count() > 0):
                if (pd.notnull(train_station['Stop:Arrival time'])):
                    time_diffs = abs(nearest_weather_station['start_time'] - train_station['Stop:Arrival time'])

                    nearest_time_index = time_diffs.idxmin()
                    nearest_weather_data = nearest_weather_station.loc[nearest_time_index]

                    combined_row = pd.concat([train_station, nearest_weather_data])
                    combined_data.append(combined_row)

    combined_df = pd.DataFrame(combined_data)

    combined_df.reset_index(drop=True, inplace=True)
    
    combined_df.to_csv(CSV_save_location_full_data_weather_and_train, index=False)
else:
    combined_df = pd.read_csv(CSV_save_location_full_data_weather_and_train, low_memory=False)


display(combined_df)
display(combined_df.dtypes)

	Service:RDT-ID	Service:Date	Service:Type	Service:Company	Service:Train number	Service:Completely cancelled	Service:Partly cancelled	Service:Maximum delay	Stop:RDT-ID	Stop:Station code	...	VV	N	U	WW	IX	M	R	S	O	Y
0	10368001	2023-01-01	sprinter	ns	4026	False	False	0	93070750	AC	...	73.0	8.0	74.0	NaN	5	0.0	0.0	0.0	0.0	0.0
1	10368079	2023-01-01	sprinter	ns	4037	False	False	0	93071464	AC	...	74.0	8.0	74.0	NaN	5	0.0	0.0	0.0	0.0	0.0
2	10368089	2023-01-01	sprinter	ns	4035	False	False	0	93071556	AC	...	73.0	8.0	74.0	NaN	5	0.0	0.0	0.0	0.0	0.0
3	10368125	2023-01-01	sprinter	ns	4030	False	False	0	93071943	AC	...	74.0	8.0	74.0	NaN	5	0.0	0.0	0.0	0.0	0.0
4	10368142	2023-01-01	sprinter	ns	4028	False	False	0	93072094	AC	...	73.0	8.0	74.0	NaN	5	0.0	0.0	0.0	0.0	0.0
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
4288666	10894127	2023-03-28	stoptrein	blauwnet	8068	False	False	0	98003158	NaN	...	84.0	4.0	58.0	NaN	5	0.0	0.0	0.0	0.0	0.0
4288667	10894200	2023-03-28	stoptrein	blauwnet	8055	False	False	0	98003859	NaN	...	83.0	8.0	64.0	51.0	7	0.0	1.0	0.0	0.0	0.0
4288668	10894311	2023-03-28	stoptrein	blauwnet	8070	False	False	0	98004967	NaN	...	83.0	8.0	64.0	51.0	7	0.0	1.0	0.0	0.0	0.0
4288669	10894382	2023-03-28	stoptrein	blauwnet	8057	False	False	0	98005654	NaN	...	80.0	8.0	69.0	81.0	7	0.0	1.0	0.0	0.0	0.0
4288670	10894488	2023-03-28	stoptrein	blauwnet	8072	False	False	0	98006751	NaN	...	83.0	8.0	64.0	51.0	7	0.0	1.0	0.0	0.0	0.0

4288671 rows × 54 columns

Service:RDT-ID                    int64
Service:Date                     object
Service:Type                     object
Service:Company                  object
Service:Train number              int64
Service:Completely cancelled       bool
Service:Partly cancelled           bool
Service:Maximum delay             int64
Stop:RDT-ID                       int64
Stop:Station code                object
Stop:Station name                object
Stop:Arrival time                object
Stop:Arrival delay              float64
Stop:Arrival cancelled             bool
Stop:Departure time              object
Stop:Departure delay            float64
Stop:Departure cancelled         object
id                              float64
code                             object
uic                             float64
name_short                       object
name_medium                      object
name_long                        object
slug                             object
country                          object
type                             object
geo_lat                         float64
geo_lng                         float64
geo_lat_radians                 float64
geo_lng_radians                 float64
STN                               int64
start_time                       object
DD                              float64
FH                              float64
FF                              float64
FX                              float64
temp                            float64
T10N                            float64
TD                              float64
SQ                              float64
Q                               float64
DR                              float64
RH                              float64
P                               float64
VV                              float64
N                               float64
U                               float64
WW                              float64
IX                                int64
M                               float64
R                               float64
S                               float64
O                               float64
Y                               float64
dtype: object

even de tijden naar een datetime object zetten.

combined_df['Service:Date'] = pd.to_datetime(combined_df['Service:Date']).dt.dayofweek
combined_df['Stop:Departure time'] = pd.to_datetime(combined_df['Stop:Departure time'], utc=True).dt.dayofweek
combined_df['Stop:Arrival time'] = pd.to_datetime(combined_df['Stop:Arrival time'], utc=True).dt.dayofweek

Hieronder verwijder ik de rijen waar de het type van het verkeer niet een goede verandering heeft of een infinate variable in hun tests.
Deze testen zijn later in dit document te vinden, in het copje in en out control

combined_df = combined_df[~combined_df['Service:Type'].isin(['stopbus i.p.v. trein', 'snelbus i.p.v. trein', 'int. trein', 'bus'])]

Hieronder ga ik nog een paar null variables een betere betekenis geven. dit is om er voor te zorgen dat ik een beteren datapatroon krijg

combined_df = combined_df[combined_df['id'].notna()]
combined_df['Stop:Station code'] = combined_df['Stop:Station code'].fillna('NaN')
combined_df['code'] = combined_df['code'].fillna('NaN')
display(combined_df.isna().sum())

Service:RDT-ID                        0
Service:Date                          0
Service:Type                          0
Service:Company                       0
Service:Train number                  0
Service:Completely cancelled          0
Service:Partly cancelled              0
Service:Maximum delay                 0
Stop:RDT-ID                           0
Stop:Station code                     0
Stop:Station name                     0
Stop:Arrival time                     0
Stop:Arrival delay                    0
Stop:Arrival cancelled                0
Stop:Departure time              476147
Stop:Departure delay             476147
Stop:Departure cancelled         476147
id                                    0
code                                  0
uic                                   0
name_short                            0
name_medium                           0
name_long                             0
slug                                  0
country                               0
type                                  0
geo_lat                               0
geo_lng                               0
geo_lat_radians                       0
geo_lng_radians                       0
STN                                   0
start_time                            0
DD                                92827
FH                                92827
FF                                92827
FX                                92827
temp                             149644
T10N                            3538233
TD                               149644
SQ                               230628
Q                                230628
DR                               149722
RH                               149644
P                                963340
VV                               875951
N                               1168562
U                                149644
WW                              2653476
IX                                    0
M                                879361
R                                879361
S                                879361
O                                879361
Y                                879361
dtype: int64

Hieronder ga ik controleren of de data overeenkomt met wat ik verwacht.
Dit betekent dat de geografische coördinaten en de tijden dicht bij elkaar moeten liggen.

combined_df = combined_df.merge(weather_station_location, on='STN', how='inner')
display(combined_df[['LAT(north)', 'LON(east)','geo_lat', 'geo_lng', 'start_time', 'Stop:Arrival time']])
display(combined_df.info())

	LAT(north)	LON(east)	geo_lat	geo_lng	start_time	Stop:Arrival time
0	52.318	4.790	52.278500	4.977000	2023-01-01 10:00:00+01:00	6
1	52.318	4.790	52.278500	4.977000	2023-01-01 11:00:00+01:00	6
2	52.318	4.790	52.278500	4.977000	2023-01-01 10:00:00+01:00	6
3	52.318	4.790	52.278500	4.977000	2023-01-01 11:00:00+01:00	6
4	52.318	4.790	52.278500	4.977000	2023-01-01 10:00:00+01:00	6
...	...	...	...	...	...	...
4179565	52.750	6.574	52.718613	6.848611	2023-03-28 17:00:00+02:00	1
4179566	52.750	6.574	52.718613	6.848611	2023-03-28 18:00:00+02:00	1
4179567	52.750	6.574	52.718613	6.848611	2023-03-28 18:00:00+02:00	1
4179568	52.750	6.574	52.718613	6.848611	2023-03-28 19:00:00+02:00	1
4179569	52.750	6.574	52.718613	6.848611	2023-03-28 18:00:00+02:00	1

4179570 rows × 6 columns

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 4179570 entries, 0 to 4179569
Data columns (total 58 columns):
 #   Column                        Dtype  
---  ------                        -----  
 0   Service:RDT-ID                int64  
 1   Service:Date                  int32  
 2   Service:Type                  object 
 3   Service:Company               object 
 4   Service:Train number          int64  
 5   Service:Completely cancelled  bool   
 6   Service:Partly cancelled      bool   
 7   Service:Maximum delay         int64  
 8   Stop:RDT-ID                   int64  
 9   Stop:Station code             object 
 10  Stop:Station name             object 
 11  Stop:Arrival time             int32  
 12  Stop:Arrival delay            float64
 13  Stop:Arrival cancelled        bool   
 14  Stop:Departure time           float64
 15  Stop:Departure delay          float64
 16  Stop:Departure cancelled      object 
 17  id                            float64
 18  code                          object 
 19  uic                           float64
 20  name_short                    object 
 21  name_medium                   object 
 22  name_long                     object 
 23  slug                          object 
 24  country                       object 
 25  type                          object 
 26  geo_lat                       float64
 27  geo_lng                       float64
 28  geo_lat_radians               float64
 29  geo_lng_radians               float64
 30  STN                           int64  
 31  start_time                    object 
 32  DD                            float64
 33  FH                            float64
 34  FF                            float64
 35  FX                            float64
 36  temp                          float64
 37  T10N                          float64
 38  TD                            float64
 39  SQ                            float64
 40  Q                             float64
 41  DR                            float64
 42  RH                            float64
 43  P                             float64
 44  VV                            float64
 45  N                             float64
 46  U                             float64
 47  WW                            float64
 48  IX                            int64  
 49  M                             float64
 50  R                             float64
 51  S                             float64
 52  O                             float64
 53  Y                             float64
 54  LON(east)                     float64
 55  LAT(north)                    float64
 56  ALT(m)                        float64
 57  NAME                          object 
dtypes: bool(3), float64(33), int32(2), int64(6), object(14)
memory usage: 1.7+ GB

None

Uit de gegevens hierboven kan ik afleiden dat het gelukt is om het dichtstbijzijnde weerstation te vinden voor een treinstation, met de tijd als extra parameter.
Met deze data kan ik eindelijk beginnen met het analyseren van wat een machine learning algoritme denkt over de data. Maar eerst ga ik kijken hoe de correlaties van deze hele dataset eruitzien.
Voordat ik dat doe, moet ik echter enkele kolommen verwijderen die dubbel zijn of niet nodig.
En dat is precies wat ik hieronder ga doen.

combined_df = combined_df.drop(columns=['name_medium', 'name_short', 'name_long', 'LAT(north)', 'LON(east)','geo_lat', 'geo_lng', 'NAME', 'ALT(m)', 'start_time' , 'STN', 'geo_lat_radians', 'geo_lng_radians',
                                                   'code', 'Stop:Station code', 'id',  'Service:RDT-ID', 'Service:Train number', 'Service:Company', 'slug'])
display(combined_df.info())

<class 'pandas.core.frame.DataFrame'>
RangeIndex: 4179570 entries, 0 to 4179569
Data columns (total 38 columns):
 #   Column                        Dtype  
---  ------                        -----  
 0   Service:Date                  int32  
 1   Service:Type                  object 
 2   Service:Completely cancelled  bool   
 3   Service:Partly cancelled      bool   
 4   Service:Maximum delay         int64  
 5   Stop:RDT-ID                   int64  
 6   Stop:Station name             object 
 7   Stop:Arrival time             int32  
 8   Stop:Arrival delay            float64
 9   Stop:Arrival cancelled        bool   
 10  Stop:Departure time           float64
 11  Stop:Departure delay          float64
 12  Stop:Departure cancelled      object 
 13  uic                           float64
 14  country                       object 
 15  type                          object 
 16  DD                            float64
 17  FH                            float64
 18  FF                            float64
 19  FX                            float64
 20  temp                          float64
 21  T10N                          float64
 22  TD                            float64
 23  SQ                            float64
 24  Q                             float64
 25  DR                            float64
 26  RH                            float64
 27  P                             float64
 28  VV                            float64
 29  N                             float64
 30  U                             float64
 31  WW                            float64
 32  IX                            int64  
 33  M                             float64
 34  R                             float64
 35  S                             float64
 36  O                             float64
 37  Y                             float64
dtypes: bool(3), float64(25), int32(2), int64(3), object(5)
memory usage: 1.1+ GB

None

Hieronder maak ik de correlatieplot die hopelijk zal helpen om de correlaties te ontdekken.

plt.subplots(figsize=(15, 15))
corr = combined_df.select_dtypes(include=['int32', 'int64', 'float64']).corr()
mask = np.tril(corr)
sns.heatmap(corr, mask=mask, cbar=True, annot=True, fmt=".1f")

<Axes: >

Uit de correlatieplot is al redelijk wat zichtbaar. Maar om het duidelijker te maken, ga ik ze allemaal met elkaar vergelijken in plots.

Hieronder creëer ik de plots om dit goed te laten zien.

plt.figure(figsize=(40,40))

# sns.pairplot(combined_df.sample(frac=0.3)[['Service:Type','Stop:Arrival time', 'Stop:Departure time', 'Stop:Departure delay', 'Stop:Arrival delay', 'temp', 'TD', 'SQ', 'Q', 'DR', 'RH', 'P', 'VV', 'N', 'U', 'IX', 'M', 'R', 'S', 'O', 'Y']], hue='Service:Type', corner= True, palette='rainbow')
input = combined_df[combined_df['country'] == 'NL']

Yfeatures = ['Stop:Departure delay', 'Stop:Arrival delay', 'Service:Maximum delay', 'Stop:Arrival cancelled', 'Stop:Departure cancelled', 'Service:Partly cancelled', 'Service:Completely cancelled'] # combined_df.columns
Xfeatures = ['DD' ,'FH' ,'FF' ,'FX' ,'temp' ,'T10N' ,'TD' ,'SQ' ,'Q' ,'DR' ,'RH' ,'P' ,'VV' ,'N' ,'U' ,'WW' ,'M' ,'R' ,'S' ,'O' ,'Y'] # combined_df.columns

pcount=1

plt.figure(figsize=(60, 20))
for Y in Yfeatures:
    for X in Xfeatures:
        plot = plt.subplot(Yfeatures.__len__(), Xfeatures.__len__(), pcount)
        sns.histplot(legend=True, color=[1,0,0], data=input, x=X, y=Y, bins=50)
        pcount +=1

plt.savefig('Combimed_Weather_Train_Pairplot.png')

<Figure size 4000x4000 with 0 Axes>

Hierboven toon ik de correlaties op een geplotte manier.
Het lijkt er echter op dat ik nog wat verder moet zoeken om te bepalen of er daadwerkelijk een manier is om dit te voorspellen.
Op dit moment lijkt het erop dat ik weinig heb om mee te werken.

📃🪓 Splitting naar train/test

het is nu tijd om mijn data te splitten tussen training data en test data.
Dit doe ik voor de testen en trainen van de machine learing modellen.
Hier gebruik ik train_test_split van sklearn voor, dit is mij aangeleerd in eerdere opdrachten.
Hier gebruik ik ook sklearn's preprocessing om ervoor te zorgen dat alle lebbels goed mee werkt.

from sklearn.model_selection import train_test_split
from pandas.api.types import is_string_dtype
from pandas.api.types import is_numeric_dtype
from pandas.api.types import is_bool_dtype
from sklearn import preprocessing


Smaller_combined_df = combined_df.sample(frac=0.3).drop(columns=['Stop:Departure delay', 'Stop:Departure cancelled'])

le = preprocessing.LabelEncoder()
scaler = preprocessing.MinMaxScaler()

for column in Smaller_combined_df.columns:
    if(is_string_dtype(Smaller_combined_df[column])): 
        Smaller_combined_df[column] = le.fit_transform(Smaller_combined_df[column])
    else: 
        if (is_numeric_dtype(Smaller_combined_df[column]) and not is_bool_dtype(Smaller_combined_df[column])):
            Smaller_combined_df[column] = scaler.fit_transform(Smaller_combined_df[column].values.reshape(-1, 1))

Smaller_combined_df = Smaller_combined_df.fillna(-1)

X = Smaller_combined_df.drop(columns=['Stop:Arrival delay', 'Service:Maximum delay'])
y = Smaller_combined_df[['Stop:Arrival delay']]

X_train, X_test, y_train, y_test = train_test_split( X, y, random_state=21 )
print("There are in total", len(X), "observations, of which", len(X_train), "are now in the train set, and", len(X_test), "in the test set.")
display(X_train.dtypes)
display(X_train.isna().sum())
display(X_train.info())

There are in total 1253871 observations, of which 940403 are now in the train set, and 313468 in the test set.

Service:Date                    float64
Service:Type                      int32
Service:Completely cancelled       bool
Service:Partly cancelled           bool
Stop:RDT-ID                     float64
Stop:Station name                 int32
Stop:Arrival time               float64
Stop:Arrival cancelled             bool
Stop:Departure time             float64
uic                             float64
country                           int32
type                              int32
DD                              float64
FH                              float64
FF                              float64
FX                              float64
temp                            float64
T10N                            float64
TD                              float64
SQ                              float64
Q                               float64
DR                              float64
RH                              float64
P                               float64
VV                              float64
N                               float64
U                               float64
WW                              float64
IX                              float64
M                               float64
R                               float64
S                               float64
O                               float64
Y                               float64
dtype: object

Service:Date                    0
Service:Type                    0
Service:Completely cancelled    0
Service:Partly cancelled        0
Stop:RDT-ID                     0
Stop:Station name               0
Stop:Arrival time               0
Stop:Arrival cancelled          0
Stop:Departure time             0
uic                             0
country                         0
type                            0
DD                              0
FH                              0
FF                              0
FX                              0
temp                            0
T10N                            0
TD                              0
SQ                              0
Q                               0
DR                              0
RH                              0
P                               0
VV                              0
N                               0
U                               0
WW                              0
IX                              0
M                               0
R                               0
S                               0
O                               0
Y                               0
dtype: int64

<class 'pandas.core.frame.DataFrame'>
Index: 940403 entries, 794871 to 3519301
Data columns (total 34 columns):
 #   Column                        Non-Null Count   Dtype  
---  ------                        --------------   -----  
 0   Service:Date                  940403 non-null  float64
 1   Service:Type                  940403 non-null  int32  
 2   Service:Completely cancelled  940403 non-null  bool   
 3   Service:Partly cancelled      940403 non-null  bool   
 4   Stop:RDT-ID                   940403 non-null  float64
 5   Stop:Station name             940403 non-null  int32  
 6   Stop:Arrival time             940403 non-null  float64
 7   Stop:Arrival cancelled        940403 non-null  bool   
 8   Stop:Departure time           940403 non-null  float64
 9   uic                           940403 non-null  float64
 10  country                       940403 non-null  int32  
 11  type                          940403 non-null  int32  
 12  DD                            940403 non-null  float64
 13  FH                            940403 non-null  float64
 14  FF                            940403 non-null  float64
 15  FX                            940403 non-null  float64
 16  temp                          940403 non-null  float64
 17  T10N                          940403 non-null  float64
 18  TD                            940403 non-null  float64
 19  SQ                            940403 non-null  float64
 20  Q                             940403 non-null  float64
 21  DR                            940403 non-null  float64
 22  RH                            940403 non-null  float64
 23  P                             940403 non-null  float64
 24  VV                            940403 non-null  float64
 25  N                             940403 non-null  float64
 26  U                             940403 non-null  float64
 27  WW                            940403 non-null  float64
 28  IX                            940403 non-null  float64
 29  M                             940403 non-null  float64
 30  R                             940403 non-null  float64
 31  S                             940403 non-null  float64
 32  O                             940403 non-null  float64
 33  Y                             940403 non-null  float64
dtypes: bool(3), float64(27), int32(4)
memory usage: 217.9 MB

None

het ziet er naar uit dat de data hierboven goed is gesplit.
Er is alleen wijnig in te zien waar een draad an te knowpen is.
Maar nog steeds zal ik kijken wat er in de sets zit hieronder.

display(X_train.tail(5))
display(y_train.tail(5))

	Service:Date	Service:Type	Service:Completely cancelled	Service:Partly cancelled	Stop:RDT-ID	Stop:Station name	Stop:Arrival time	Stop:Arrival cancelled	Stop:Departure time	uic	...	VV	N	U	WW	IX	M	R	S	O	Y
3526960	0.000000	3	False	False	0.841239	446	0.000000	False	-1.000000	0.756620	...	0.285714	0.888889	0.936709	0.565217	1.0	0.0	1.0	0.0	0.0	0.0
1239497	1.000000	9	False	False	0.754576	135	1.000000	False	1.000000	0.756382	...	-1.000000	-1.000000	0.708861	-1.000000	0.5	-1.0	-1.0	-1.0	-1.0	-1.0
164381	0.666667	9	False	False	0.885239	20	0.666667	False	0.666667	0.756338	...	0.952381	0.888889	0.569620	0.250000	1.0	0.0	1.0	0.0	0.0	0.0
4012454	0.666667	3	False	False	0.431312	510	0.666667	False	0.666667	0.756680	...	0.928571	0.888889	0.607595	-1.000000	0.0	0.0	0.0	0.0	0.0	0.0
564179	0.833333	9	False	False	0.897767	32	0.833333	False	0.833333	0.756313	...	0.833333	0.888889	0.772152	0.250000	1.0	0.0	1.0	0.0	0.0	0.0

5 rows × 34 columns

	Stop:Arrival delay
3526960	0.000000
1239497	0.000000
164381	0.008929
4012454	0.000000
564179	0.000000

display(X_test.head(5))
display(y_test.head(5))

	Service:Date	Service:Type	Service:Completely cancelled	Service:Partly cancelled	Stop:RDT-ID	Stop:Station name	Stop:Arrival time	Stop:Arrival cancelled	Stop:Departure time	uic	...	VV	N	U	WW	IX	M	R	S	O	Y
1319171	1.000000	10	False	False	0.297529	81	1.000000	False	1.000000	0.549328	...	-1.000000	-1.000000	0.721519	-1.00	0.5	-1.0	-1.0	-1.0	-1.0	-1.0
1794098	0.166667	3	False	False	0.246390	126	0.166667	False	0.166667	0.756434	...	0.726190	0.888889	0.810127	-1.00	0.0	0.0	0.0	0.0	0.0	0.0
2053586	0.500000	3	False	False	0.043974	208	0.500000	False	0.500000	0.756437	...	-1.000000	-1.000000	-1.000000	-1.00	0.5	-1.0	-1.0	-1.0	-1.0	-1.0
2878968	0.333333	9	False	False	0.946618	358	0.333333	False	0.333333	0.756548	...	0.952381	0.777778	0.594937	-1.00	0.0	0.0	0.0	0.0	0.0	0.0
1005833	1.000000	10	False	False	0.753972	120	1.000000	False	1.000000	0.756373	...	0.738095	0.888889	0.911392	0.25	1.0	0.0	1.0	0.0	0.0	0.0

5 rows × 34 columns

	Stop:Arrival delay
1319171	0.000000
1794098	0.000000
2053586	0.000000
2878968	0.004464
1005833	0.000000

Het ziet er uit dat het goed is gesplit en dat ik vereder kan naar het maken van een eerste poging in Machine learining.

🧬 Modelling

Decision Tree

Hieronder ga ik beginnen met een decisiontreeregressor van sklearn, dit is om te zien wat er gebeurt.

from sklearn.tree import DecisionTreeRegressor
from sklearn.metrics import mean_absolute_error,max_error,mean_squared_error,mean_squared_log_error,r2_score

tree = DecisionTreeRegressor(random_state=25, max_depth=5)

tree.fit(X_train, y_train)
tree.score(X_test, y_test)

pred = tree.predict(X_test)

msq = mean_squared_error(y_test,pred)
print("mean squared err:", msq)
mab = mean_absolute_error(y_test,pred)
print("mean absolute err:", mab)
msql = mean_squared_log_error(y_test,pred)
print("mean squared log err:", msql)
mx = max_error(y_test,pred)
print("max err:", mx)
Rsq = r2_score(y_test,pred)
print("R squared:", Rsq)

plt.scatter(pred, y_test)
plt.xlabel('predicted')
plt.ylabel('eigenlijke variable')

mean squared err: 0.00015172037582403077
mean absolute err: 0.004999627870634173
mean squared log err: 0.00012859279405038072
max err: 0.9426387461459403
R squared: 0.035583038496755925

Text(0, 0.5, 'eigenlijke variable')

Het ziet ernaaruit dat hierboven er geen goede verbinding is om te voorspellen hoelang een vertraging is te bepalen, Want de test gaat er vanuit dat er niet zo veel vertraging is.
zelfs all is er erge weeromstandigheden gemeld.

hieronder laat ik zien waarom ik hyper parameters gebruik die ik gebruik hierboven.

from sklearn.model_selection import GridSearchCV
parameters = {'max_depth':list(range(5, 21)), 'random_state':list(range(20, 30))}

clf = GridSearchCV(DecisionTreeRegressor(), parameters, n_jobs=10)
clf.fit(X_train, y_train)

df = pd.DataFrame(clf.cv_results_)
df.sort_values(by='rank_test_score')

	mean_fit_time	std_fit_time	mean_score_time	std_score_time	param_max_depth	param_random_state	params	split0_test_score	split1_test_score	split2_test_score	split3_test_score	split4_test_score	mean_test_score	std_test_score	rank_test_score
9	1.383550	0.052501	0.032910	0.004670	5	29	{'max_depth': 5, 'random_state': 29}	0.032087	0.036419	0.076407	0.040098	0.019564	0.040915	0.019049	1
8	1.360434	0.022852	0.022548	0.003704	5	28	{'max_depth': 5, 'random_state': 28}	0.032087	0.036419	0.076407	0.040098	0.019564	0.040915	0.019049	1
7	1.319590	0.022095	0.030718	0.004214	5	27	{'max_depth': 5, 'random_state': 27}	0.032087	0.036419	0.076407	0.040098	0.019564	0.040915	0.019049	3
2	1.394482	0.025924	0.023387	0.004725	5	22	{'max_depth': 5, 'random_state': 22}	0.032087	0.036419	0.076407	0.040098	0.019564	0.040915	0.019049	4
5	1.320968	0.017483	0.029929	0.004009	5	25	{'max_depth': 5, 'random_state': 25}	0.032087	0.036419	0.076407	0.040098	0.019564	0.040915	0.019049	4
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
150	4.040698	0.030783	0.027346	0.003383	20	20	{'max_depth': 20, 'random_state': 20}	-0.238856	-0.407118	-0.488603	-0.336659	-0.527773	-0.399802	0.104122	156
151	4.084840	0.032983	0.027395	0.001086	20	21	{'max_depth': 20, 'random_state': 21}	-0.280658	-0.540218	-0.298487	-0.316655	-0.569125	-0.401028	0.126296	157
159	3.915019	0.037795	0.021045	0.001988	20	29	{'max_depth': 20, 'random_state': 29}	-0.258677	-0.435636	-0.467134	-0.315934	-0.599047	-0.415285	0.119388	158
154	4.021103	0.018816	0.026219	0.003296	20	24	{'max_depth': 20, 'random_state': 24}	-0.291595	-0.479860	-0.445227	-0.345751	-0.520413	-0.416569	0.085140	159
157	4.023842	0.025901	0.027158	0.004099	20	27	{'max_depth': 20, 'random_state': 27}	-0.238492	-0.522900	-0.439356	-0.395254	-0.523385	-0.423878	0.105008	160

160 rows × 15 columns

hieronder zal ik de decisionTree weergeven zodat het zichtbaar is wat het doet of meer gezien wat het waarschijnlijk doet.

from sklearn.tree import plot_tree
plt.figure(figsize=(180,100))
plot_tree(tree, fontsize=8, feature_names=X.columns)

[Text(0.5823170731707317, 0.9166666666666666, 'Service:Type <= 13.5\nsquared_error = 0.0\nsamples = 321632\nvalue = 0.004'),
 Text(0.3048780487804878, 0.75, 'Service:Type <= 3.5\nsquared_error = 0.0\nsamples = 320799\nvalue = 0.004'),
 Text(0.13414634146341464, 0.5833333333333334, 'Stop:Station name <= 26.5\nsquared_error = 0.005\nsamples = 1019\nvalue = 0.026'),
 Text(0.04878048780487805, 0.4166666666666667, 'DD <= 0.03\nsquared_error = 0.016\nsamples = 99\nvalue = 0.073'),
 Text(0.024390243902439025, 0.25, 'squared_error = 0.0\nsamples = 1\nvalue = 0.519'),
 Text(0.07317073170731707, 0.25, 'DD <= 0.217\nsquared_error = 0.015\nsamples = 98\nvalue = 0.069'),
 Text(0.04878048780487805, 0.08333333333333333, 'squared_error = 0.024\nsamples = 48\nvalue = 0.11'),
 Text(0.0975609756097561, 0.08333333333333333, 'squared_error = 0.003\nsamples = 50\nvalue = 0.029'),
 Text(0.21951219512195122, 0.4166666666666667, 'Service:Maximum delay <= 0.1\nsquared_error = 0.003\nsamples = 920\nvalue = 0.021'),
 Text(0.17073170731707318, 0.25, 'RH <= 0.418\nsquared_error = 0.003\nsamples = 917\nvalue = 0.021'),
 Text(0.14634146341463414, 0.08333333333333333, 'squared_error = 0.003\nsamples = 915\nvalue = 0.02'),
 Text(0.1951219512195122, 0.08333333333333333, 'squared_error = 0.041\nsamples = 2\nvalue = 0.204'),
 Text(0.2682926829268293, 0.25, 'Service:Date <= 0.417\nsquared_error = 0.02\nsamples = 3\nvalue = 0.201'),
 Text(0.24390243902439024, 0.08333333333333333, 'squared_error = 0.001\nsamples = 2\nvalue = 0.301'),
 Text(0.2926829268292683, 0.08333333333333333, 'squared_error = 0.0\nsamples = 1\nvalue = 0.0'),
 Text(0.47560975609756095, 0.5833333333333334, 'Stop:Arrival cancelled <= 0.5\nsquared_error = 0.0\nsamples = 319780\nvalue = 0.004'),
 Text(0.4146341463414634, 0.4166666666666667, 'Service:Partly cancelled <= 0.5\nsquared_error = 0.0\nsamples = 293728\nvalue = 0.004'),
 Text(0.36585365853658536, 0.25, 'Service:Type <= 9.5\nsquared_error = 0.0\nsamples = 281197\nvalue = 0.004'),
 Text(0.34146341463414637, 0.08333333333333333, 'squared_error = 0.0\nsamples = 71975\nvalue = 0.005'),
 Text(0.3902439024390244, 0.08333333333333333, 'squared_error = 0.0\nsamples = 209222\nvalue = 0.003'),
 Text(0.4634146341463415, 0.25, 'Stop:Station name <= 404.5\nsquared_error = 0.001\nsamples = 12531\nvalue = 0.009'),
 Text(0.43902439024390244, 0.08333333333333333, 'squared_error = 0.001\nsamples = 9080\nvalue = 0.011'),
 Text(0.4878048780487805, 0.08333333333333333, 'squared_error = 0.001\nsamples = 3451\nvalue = 0.007'),
 Text(0.5365853658536586, 0.4166666666666667, 'U <= 0.943\nsquared_error = 0.0\nsamples = 26052\nvalue = 0.0'),
 Text(0.5121951219512195, 0.25, 'squared_error = 0.0\nsamples = 22871\nvalue = 0.0'),
 Text(0.5609756097560976, 0.25, 'VV <= 0.768\nsquared_error = 0.0\nsamples = 3181\nvalue = 0.0'),
 Text(0.5365853658536586, 0.08333333333333333, 'squared_error = 0.0\nsamples = 3043\nvalue = 0.0'),
 Text(0.5853658536585366, 0.08333333333333333, 'squared_error = 0.0\nsamples = 138\nvalue = 0.0'),
 Text(0.8597560975609756, 0.75, 'Service:Maximum delay <= 0.1\nsquared_error = 0.005\nsamples = 833\nvalue = 0.03'),
 Text(0.7682926829268293, 0.5833333333333334, 'uic <= 0.757\nsquared_error = 0.005\nsamples = 830\nvalue = 0.029'),
 Text(0.7073170731707317, 0.4166666666666667, 'VV <= 0.351\nsquared_error = 0.009\nsamples = 258\nvalue = 0.044'),
 Text(0.6585365853658537, 0.25, 'VV <= 0.333\nsquared_error = 0.063\nsamples = 11\nvalue = 0.16'),
 Text(0.6341463414634146, 0.08333333333333333, 'squared_error = 0.007\nsamples = 10\nvalue = 0.085'),
 Text(0.6829268292682927, 0.08333333333333333, 'squared_error = 0.0\nsamples = 1\nvalue = 0.912'),
 Text(0.7560975609756098, 0.25, 'temp <= 0.818\nsquared_error = 0.006\nsamples = 247\nvalue = 0.039'),
 Text(0.7317073170731707, 0.08333333333333333, 'squared_error = 0.006\nsamples = 242\nvalue = 0.036'),
 Text(0.7804878048780488, 0.08333333333333333, 'squared_error = 0.021\nsamples = 5\nvalue = 0.194'),
 Text(0.8292682926829268, 0.4166666666666667, 'VV <= 0.018\nsquared_error = 0.002\nsamples = 572\nvalue = 0.022'),
 Text(0.8048780487804879, 0.25, 'squared_error = 0.0\nsamples = 1\nvalue = 0.292'),
 Text(0.8536585365853658, 0.25, 'N <= -0.5\nsquared_error = 0.002\nsamples = 571\nvalue = 0.021'),
 Text(0.8292682926829268, 0.08333333333333333, 'squared_error = 0.002\nsamples = 226\nvalue = 0.013'),
 Text(0.8780487804878049, 0.08333333333333333, 'squared_error = 0.002\nsamples = 345\nvalue = 0.027'),
 Text(0.9512195121951219, 0.5833333333333334, 'VV <= 0.893\nsquared_error = 0.005\nsamples = 3\nvalue = 0.273'),
 Text(0.926829268292683, 0.4166666666666667, 'WW <= -0.223\nsquared_error = 0.001\nsamples = 2\nvalue = 0.225'),
 Text(0.9024390243902439, 0.25, 'squared_error = 0.0\nsamples = 1\nvalue = 0.199'),
 Text(0.9512195121951219, 0.25, 'squared_error = -0.0\nsamples = 1\nvalue = 0.25'),
 Text(0.975609756097561, 0.4166666666666667, 'squared_error = 0.0\nsamples = 1\nvalue = 0.37')]

hieronder laat ik de shaply van mij decision tree zien met behulp van de shap library. dit doe ik om te kijken wat er gebruikt word.

import shap

explainer = shap.Explainer(tree)
shap_values = explainer.shap_values(X_test)
# plt.figure(figsize=(18000,10000))

# shap.bar_plot(shap_values, X_test)
shap.summary_plot(shap_values, X_test)

C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\tqdm\auto.py:21: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html
  from .autonotebook import tqdm as notebook_tqdm

Hieronder ga ik kijken naar een verzie van gridsearch om te kijken welke hyperparameterers goed zijn voor de modelen die ik later ga gebruiken om een beteren accarucy eruit krijgen.

Hieronder doe ik de Random forest regressor met behulp van de Scikit-learn library hun BayesSearchCV.
het zou met deze sneler zijn om je hyperparameters te vinden volgens een andere student.

Als eerst ga ik mijn RandomForestRegressor doen.

from skopt import BayesSearchCV
from sklearn.ensemble import RandomForestRegressor

searchSpaces = { 
    'max_depth':list(range(5, 21)),
    'n_estimators':list(range(5, 16)),
    'random_state':list(range(15, 26)),
    }

opt = BayesSearchCV(RandomForestRegressor(), search_spaces=searchSpaces, n_jobs=-1, n_points=24)

opt.fit(X_train, y_train)

print("val. score: %s" % opt.best_score_)
print("test score: %s" % opt.score(X_test, y_test))
print("best params: %s" % str(opt.best_params_))

val. score: 0.035911527237130335
test score: -0.42842673041392776
best params: OrderedDict([('max_depth', 8), ('n_estimators', 15), ('random_state', 24)])

C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\base.py:1474: DataConversionWarning: A column-vector y was passed when a 1d array was expected. Please change the shape of y to (n_samples,), for example using ravel().
  return fit_method(estimator, *args, **kwargs)

Nu ga ik mijn Lasso doen.

from skopt import BayesSearchCV
from sklearn.linear_model import Lasso

searchSpaces = { 
    'alpha':np.arange(0.2, 5.0, 0.2),
    'max_iter':np.arange(100, 1500, 100),
    'random_state':list(range(15, 26)),
    }

opt = BayesSearchCV(Lasso(), search_spaces=searchSpaces, n_jobs=-1, n_points=24)

opt.fit(X_train, y_train.squeeze())

print("val. score: %s" % opt.best_score_)
print("test score: %s" % opt.score(X_test, y_test.squeeze()))
print("best params: %s" % str(opt.best_params_))

val. score: -0.00035203268165866496
test score: 0.00029357364734106106
best params: OrderedDict([('alpha', 0.2), ('max_iter', 900), ('random_state', 20)])

Nu ga ik mijn BaggingRegressor doen.

from skopt import BayesSearchCV
from sklearn.ensemble import BaggingRegressor

searchSpaces = { 
    'n_estimators':np.arange(1, 20, 1),
    'max_samples':np.arange(0.01, 1, 0.01),
    'random_state':list(range(15, 26)),
    }

opt = BayesSearchCV(BaggingRegressor(), search_spaces=searchSpaces, n_jobs=-1, n_points=24)

opt.fit(X_train, y_train.squeeze())

print("val. score: %s" % opt.best_score_)
print("test score: %s" % opt.score(X_test, y_test.squeeze()))
print("best params: %s" % str(opt.best_params_))

val. score: 0.035199834231076865
test score: 0.07224559784689166
best params: OrderedDict([('max_samples', 0.15000000000000002), ('n_estimators', 19), ('random_state', 23)])

Nu ga ik mijn BaggingRegressor doen.

from skopt import BayesSearchCV
from sklearn.svm import LinearSVR

searchSpaces = { 
    'C':np.arange(0.1, 5, 0.1),
    'max_iter':np.arange(100, 1500, 100),
    'random_state':list(range(15, 26)),
    }

opt = BayesSearchCV(LinearSVR(), search_spaces=searchSpaces, n_jobs=-1, n_points=24)

opt.fit(X_train, y_train.squeeze())

print("val. score: %s" % opt.best_score_)
print("test score: %s" % opt.score(X_test, y_test.squeeze()))
print("best params: %s" % str(opt.best_params_))

val. score: -0.08166188625481788
test score: -0.13131061708956926
best params: OrderedDict([('C', 3.9000000000000004), ('max_iter', 500), ('random_state', 20)])

C:\Users\Anne\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\svm\_classes.py:31: FutureWarning: The default value of `dual` will change from `True` to `'auto'` in 1.5. Set the value of `dual` explicitly to suppress the warning.
  warnings.warn(
C:\Users\Anne\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\svm\_base.py:1237: ConvergenceWarning: Liblinear failed to converge, increase the number of iterations.
  warnings.warn(

Nu ga ik mijn MLPRegressor doen.

from skopt import BayesSearchCV
from sklearn.neural_network import MLPRegressor

searchSpaces = { 
    "random_state":np.arange(20, 25, 1),
    "solver":['sgd', 'adam'],
    "alpha": np.arange(0.01, 1, 0.01),
    "learning_rate": ['constant', 'invscaling', 'adaptive'],
    "max_iter": np.arange(100, 300, 20),
}

opt = BayesSearchCV(MLPRegressor(), search_spaces=searchSpaces, n_jobs=-1, n_points=24)

opt.fit(X_train, y_train.squeeze())

print("val. score: %s" % opt.best_score_)
print("test score: %s" % opt.score(X_test, y_test))
print("best params: %s" % str(opt.best_params_))

---------------------------------------------------------------------------
_RemoteTraceback                          Traceback (most recent call last)
_RemoteTraceback: 
"""
Traceback (most recent call last):
  File "C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\joblib\externals\loky\process_executor.py", line 463, in _process_worker
    r = call_item()
        ^^^^^^^^^^^
  File "C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\joblib\externals\loky\process_executor.py", line 291, in __call__
    return self.fn(*self.args, **self.kwargs)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\joblib\parallel.py", line 589, in __call__
    return [func(*args, **kwargs)
           ^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\joblib\parallel.py", line 589, in <listcomp>
    return [func(*args, **kwargs)
            ^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\utils\parallel.py", line 129, in __call__
    return self.function(*args, **kwargs)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\model_selection\_validation.py", line 895, in _fit_and_score
    estimator.fit(X_train, y_train, **fit_params)
  File "C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\base.py", line 1474, in wrapper
    return fit_method(estimator, *args, **kwargs)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\neural_network\_multilayer_perceptron.py", line 752, in fit
    return self._fit(X, y, incremental=False)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\super\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\neural_network\_multilayer_perceptron.py", line 496, in _fit
    raise ValueError(
ValueError: Solver produced non-finite parameter weights. The input data may contain large values and need to be preprocessed.
"""

The above exception was the direct cause of the following exception:

ValueError                                Traceback (most recent call last)
Cell In[142], line 14
      4 searchSpaces = { 
      5     "random_state":np.arange(20, 25, 1),
      6     "solver":['sgd', 'adam'],
   (...)
      9     "max_iter": np.arange(100, 300, 20),
     10 }
     12 opt = BayesSearchCV(MLPRegressor(), search_spaces=searchSpaces, n_jobs=-1, n_points=24)
---> 14 opt.fit(X_train, y_train)
     16 print("val. score: %s" % opt.best_score_)
     17 print("test score: %s" % opt.score(X_test, y_test))

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\skopt\searchcv.py:538, in BayesSearchCV.fit(self, X, y, groups, callback, **fit_params)
    531 if callable(self.refit):
    532     raise ValueError(
    533         "BayesSearchCV doesn't support a callable refit, "
    534         "as it doesn't define an implicit score to "
    535         "optimize"
    536     )
--> 538 super().fit(X=X, y=y, groups=groups, **fit_params)
    540 # BaseSearchCV never ranked train scores,
    541 # but apparently we used to ship this (back-compat)
    542 if self.return_train_score:

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\base.py:1474, in _fit_context.<locals>.decorator.<locals>.wrapper(estimator, *args, **kwargs)
   1467     estimator._validate_params()
   1469 with config_context(
   1470     skip_parameter_validation=(
   1471         prefer_skip_nested_validation or global_skip_validation
   1472     )
   1473 ):
-> 1474     return fit_method(estimator, *args, **kwargs)

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\model_selection\_search.py:970, in BaseSearchCV.fit(self, X, y, **params)
    964     results = self._format_results(
    965         all_candidate_params, n_splits, all_out, all_more_results
    966     )
    968     return results
--> 970 self._run_search(evaluate_candidates)
    972 # multimetric is determined here because in the case of a callable
    973 # self.scoring the return type is only known after calling
    974 first_test_score = all_out[0]["test_scores"]

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\skopt\searchcv.py:595, in BayesSearchCV._run_search(self, evaluate_candidates)
    591 while n_iter > 0:
    592     # when n_iter < n_points points left for evaluation
    593     n_points_adjusted = min(n_iter, n_points)
--> 595     optim_result, score_name = self._step(
    596         search_space,
    597         optimizer,
    598         score_name,
    599         evaluate_candidates,
    600         n_points=n_points_adjusted,
    601     )
    602     n_iter -= n_points
    604     if eval_callbacks(callbacks, optim_result):

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\skopt\searchcv.py:449, in BayesSearchCV._step(self, search_space, optimizer, score_name, evaluate_candidates, n_points)
    446 # make lists into dictionaries
    447 params_dict = [point_asdict(search_space, p) for p in params]
--> 449 all_results = evaluate_candidates(params_dict)
    451 # if self.scoring is a callable, we have to wait until here
    452 # to get the score name
    453 if score_name is None:

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\model_selection\_search.py:916, in BaseSearchCV.fit.<locals>.evaluate_candidates(candidate_params, cv, more_results)
    908 if self.verbose > 0:
    909     print(
    910         "Fitting {0} folds for each of {1} candidates,"
    911         " totalling {2} fits".format(
    912             n_splits, n_candidates, n_candidates * n_splits
    913         )
    914     )
--> 916 out = parallel(
    917     delayed(_fit_and_score)(
    918         clone(base_estimator),
    919         X,
    920         y,
    921         train=train,
    922         test=test,
    923         parameters=parameters,
    924         split_progress=(split_idx, n_splits),
    925         candidate_progress=(cand_idx, n_candidates),
    926         **fit_and_score_kwargs,
    927     )
    928     for (cand_idx, parameters), (split_idx, (train, test)) in product(
    929         enumerate(candidate_params),
    930         enumerate(cv.split(X, y, **routed_params.splitter.split)),
    931     )
    932 )
    934 if len(out) < 1:
    935     raise ValueError(
    936         "No fits were performed. "
    937         "Was the CV iterator empty? "
    938         "Were there no candidates?"
    939     )

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\sklearn\utils\parallel.py:67, in Parallel.__call__(self, iterable)
     62 config = get_config()
     63 iterable_with_config = (
     64     (_with_config(delayed_func, config), args, kwargs)
     65     for delayed_func, args, kwargs in iterable
     66 )
---> 67 return super().__call__(iterable_with_config)

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\joblib\parallel.py:1952, in Parallel.__call__(self, iterable)
   1946 # The first item from the output is blank, but it makes the interpreter
   1947 # progress until it enters the Try/Except block of the generator and
   1948 # reach the first `yield` statement. This starts the aynchronous
   1949 # dispatch of the tasks to the workers.
   1950 next(output)
-> 1952 return output if self.return_generator else list(output)

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\joblib\parallel.py:1595, in Parallel._get_outputs(self, iterator, pre_dispatch)
   1592     yield
   1594     with self._backend.retrieval_context():
-> 1595         yield from self._retrieve()
   1597 except GeneratorExit:
   1598     # The generator has been garbage collected before being fully
   1599     # consumed. This aborts the remaining tasks if possible and warn
   1600     # the user if necessary.
   1601     self._exception = True

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\joblib\parallel.py:1699, in Parallel._retrieve(self)
   1692 while self._wait_retrieval():
   1693 
   1694     # If the callback thread of a worker has signaled that its task
   1695     # triggered an exception, or if the retrieval loop has raised an
   1696     # exception (e.g. `GeneratorExit`), exit the loop and surface the
   1697     # worker traceback.
   1698     if self._aborting:
-> 1699         self._raise_error_fast()
   1700         break
   1702     # If the next job is not ready for retrieval yet, we just wait for
   1703     # async callbacks to progress.

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\joblib\parallel.py:1734, in Parallel._raise_error_fast(self)
   1730 # If this error job exists, immediatly raise the error by
   1731 # calling get_result. This job might not exists if abort has been
   1732 # called directly or if the generator is gc'ed.
   1733 if error_job is not None:
-> 1734     error_job.get_result(self.timeout)

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\joblib\parallel.py:736, in BatchCompletionCallBack.get_result(self, timeout)
    730 backend = self.parallel._backend
    732 if backend.supports_retrieve_callback:
    733     # We assume that the result has already been retrieved by the
    734     # callback thread, and is stored internally. It's just waiting to
    735     # be returned.
--> 736     return self._return_or_raise()
    738 # For other backends, the main thread needs to run the retrieval step.
    739 try:

File ~\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.11_qbz5n2kfra8p0\LocalCache\local-packages\Python311\site-packages\joblib\parallel.py:754, in BatchCompletionCallBack._return_or_raise(self)
    752 try:
    753     if self.status == TASK_ERROR:
--> 754         raise self._result
    755     return self._result
    756 finally:

ValueError: Solver produced non-finite parameter weights. The input data may contain large values and need to be preprocessed.

Nadat we hierboven alles hebben bekeken en nu weten wat de beste hyperparameters zijn voor de data.

from sklearn.neural_network import MLPRegressor
from sklearn.ensemble import BaggingRegressor
from sklearn.linear_model import Lasso
from sklearn.ensemble import RandomForestRegressor
from sklearn.linear_model import RidgeCV

from sklearn.ensemble import StackingRegressor
from sklearn.linear_model import ElasticNet


estimators = [ 
    ('mlp', MLPRegressor(alpha=0.4, learning_rate='adaptive' , max_iter=240, random_state=22, solver='adam')),
    ('brr', BaggingRegressor(max_samples=0.78, n_estimators=19, random_state=16)),
    ('las', Lasso(alpha=0.2, max_iter=1100, random_state=17)),
    ('rfr', RandomForestRegressor(n_estimators=15, random_state=22, max_depth=6)),
    ('ridge', RidgeCV()),
]

reg = StackingRegressor(estimators=estimators, final_estimator=ElasticNet())

Result = reg.fit(X_train, y_train.squeeze())

Hieronder laat ik de scores zien van de hierboven getrainde stacking model.

print(f"scoren: {reg.score(X_test, y_test.squeeze())}")

scoren: -1.1500525039620513e-05

Na de gegevens van hieboven te hebben bekeken en de nummers die ik terug krijg tot me in heb genomen.
Kan ik zeggen dat het weer niet veel te maken heeft met de hoeveelheid vertraging een trein oplopt.

Conclussie

Na het bekijken van de data die ik heb geregen en de statistieken die ik terug krijg van de testen die ik heb gedaan is het onwaarschijnlijk dat ik kan bepallen hoeveel minuten langer je zall moeten wachten op je trein.
Dit is een teleurstelling met wat ik zelf had verwacht.

in en out of controlle controllen

A/B testen

hieronder controlleer ik of de service:Type's vergenoeg uit elkaar liggen om te kijken of er iets is om mee te werken.
In de eerste code block hieronder check ik of de sprinters en de intercity's ver genoeg uit elkaar zitten om infromatie eruit te halen.

from scipy import stats

group_a = combined_df[combined_df['Service:Type'] == 'sprinter']['Stop:Arrival delay']
group_b = combined_df[combined_df['Service:Type'] == 'intercity']['Stop:Arrival delay']

t_stat, p_val = stats.ttest_ind(group_b, group_a)

print(f"T-statistic: {t_stat}")
print(f"P-value: {p_val}")

if p_val < 0.05:
    print("There is a significant difference between the two groups.")
else:
    print("There is no significant difference between the two groups.")

T-statistic: 90.92665650883148
P-value: 0.0
There is a significant difference between the two groups.

In de tweede code blok hieronder kijk ik naar all inbegrepen Service types tegen over de mean van all de Stop:Arrival delay.
Dit doe ik door alle unieke service types op te vragen en dan erdoor heen lopen, vergelijken met de mean en dan kijken of er een probleem oploopt.

import pandas as pd
from scipy import stats

service_types = combined_df['Service:Type'].unique()

results = []

for service_type in service_types:
    group = combined_df[combined_df['Service:Type'] == service_type]['Stop:Arrival delay']
    t_stat, p_val = stats.ttest_1samp(group, combined_df['Stop:Arrival delay'].mean())
    
    results.append({
        'Service:Type': service_type,
        'T-statistic': t_stat,
        'P-value': p_val,
        'Significant': 'Yes' if p_val < 0.05 else 'No'
    })

print(pd.DataFrame(results))

         Service:Type  T-statistic        P-value Significant
0            sprinter   -80.526057   0.000000e+00         Yes
1         extra trein    12.839285   7.112321e-35         Yes
2           intercity    50.374903   0.000000e+00         Yes
3              thalys    36.740517  2.111804e-279         Yes
4            nightjet    22.536041  1.980897e-108         Yes
5           stoptrein   -77.915583   0.000000e+00         Yes
6   ice international    33.593550  2.307952e-234         Yes
7           sneltrein    18.842498   4.387846e-79         Yes
8          nachttrein     5.791327   1.848929e-08         Yes
9    intercity direct    71.358031   0.000000e+00         Yes
10           eurostar    13.592323   1.910054e-40         Yes
11     speciale trein   -16.612243   1.606630e-56         Yes

Het ziet ernaar uit dat de meeste service types goed gebruikt kan worden.
Maar der zitten er een paar tussen die niet goed gebruikt kan worden en meschien weg zullen moeten gelaten worden.

Hieronder laat ik zien hoe de data eruit zag voordat ik service types verwijderd out de dataset

maar in conclusie van de A/B testen het ziet er naar uit dat met all de aanpassingen de kans op het berekenen van de tijd van vertraging is niet mogelijk of slecht mogelijk.